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Abstract 

This  investigation  explored  the  relative  performance  of  several  small-sample 
point  and  interval  estimators  for  series  system  reliability.  Among  point  estimators, 
the  maximum  likelihood  estimator  (MLE)  was  compared  to  the  corresponding  Bayes 
estimator.  In  addition,  four  interval  estimators  were  compared:  Easterling’s  mod¬ 
ified  maximum  likelihood  integer  estimator,  the  Lindstrom-Madden  estimator,  and 
Bayesian  probability  interval  estimators  constructed  using  approximate  beta  and 
Bayes  Monte  Carlo  empirical  posterior  densities. 

The  relative  performance  of  the  point  estimators  was  assessed  by  comparing 
their  mean  square  errors.  For  the  four  interval  estimators,  the  interval  coverage 
probability  and  the  average  interval  lower  bound  were  examined.  The  values  of 
these  performance  measures  were  generated  for  32  representative  series  systems  using 
Monte  Carlo  simulation. 

The  results  of  the  Monte  Carlo  study  showed  that  the  accuracy  of  the  available 
prior  information  determines  whether  a  Bayesian  or  a  classical  approach  should  be 
used  to  estimate  series  system  reliability.  If  there  is  high  confidence  that  the  mean  of 
each  component  prior  distribution  is  within  20  percent  of  the  true  component  relia¬ 
bility,  Bayesian  point  and  interval  estimators  constructed  from  an  approximate  beta 
posterior  should  be  used.  Otherwise,  the  MLE  point  estimator  and  the  Lindstrom- 
Madden  interval  estimator  should  be  chosen. 
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POINT  AND  INTERVAL  ESTIMATION  OF 
SERIES  SYSTEM  RELIABILITY  USING 
SMALL  DATA  SETS 


1.  Introduction 


Background 

The  allocation  of  firepower  to  a  list  of  key  enemy  targets  is  an  important  part 
of  strategic  war  planning.  When  assigning  scarce  strategic  resources  to  a  particular 
target,  operations  planners  need  an  accurate  assessment  of  the  likelihood  that  an 
assigned  weapon  (the  warhead  and  delivery  system)  will  neutralize  that  target.  This 
likelihood  is  called  damage  expectancy,  and  can  be  evaluated  using  the  equation 


D  =  PaXPk  (1) 

where  R,  is  the  probability  that  the  warhead  will  detonate  on  target  and  Pk  is  the 
probability  that  the  detonation  will  inflict  the  desired  damage  level.  The  quantity 
Pa  is  a  function  of  several  other  probabilities.  It  is  calculated  using  the  equation 

Pa  =  Ax  P,x  Rx  Pj,  (2) 

where  A  is  weapon  system  availability,  P,  is  the  probability  that  thf  eapon  system 
survives  launch,  R  is  the  weapon  system  reliability  (WSR)  and  Pp  is  the  probability 
that  the  weapon  system  penetrates  enemy  defenses  [2]. 

While  all  aspects  of  damage  expectancy  are  important,  this  research  focuses 
on  estimating  WSR.  According  to  Cardaronella  et  al.,  WSR  is  the  probability  that  a 


strategic  weapon  will  detonate  near  the  intended  target,  with  no  consideration  of  the 
effect  of  enemy  defenses  [7].  This  means  that,  for  a  strategic  mission  to  be  successful, 
the  delivery  system  must  get  within  an  acceptable  radius  of  the  target  center  and 
the  warhead  must  detonate  at  the  proper  yield. 

Modeling  WSR.  Headquarters  Air  Combat  Command  (ACC)  is  responsible  to 
the  Joint  Staff  for  estimating  WSR  for  weapon  systems  that  include  manned  bomber 
aircraft.  ACC’s  WSR  model  is  based  on  statistical  reliability  theory.  Since  damage 
expectancy  is  a  single-weapon  probability,  the  model  assumes  the  weapon  delivery 
system  is  composed  of  three  critical  elements:  the  bomber  aircraft,  the  air-launched 
missile,  and  a  single  warhead.  The  ACC  model  is  comprehensively  described  in  [7]; 
the  remainder  of  this  section  is  based  on  material  from  that  source  and  from  [13]. 

To  construct  the  WSR  model,  the  weapon  system  is  decomposed  into  a  series 
structure  with  three  independent  components.  These  components  correspond  to 
the  three  critical  elements  listed  above.  Each  of  these  components  can  be  viewed 
as  a  “subsystem,”  although  the  calculation  of  hardware  reliability  is  not  the  goal 
of  the  model.  Rather,  each  subsystem  in  the  reliability  model  can  be  viewed  as  the 
contribution  to  mission  success  of  each  of  the  critical  elements  of  the  weapon  system. 
The  series  model  is  used  because  if  any  of  these  elements  fails  to  complete  its  portion 
of  the  mission,  the  mission  will  not  be  a  success.  The  summary  WSR  model  is  shown 
in  Figure  1. 


Figure  1.  Summary  WSR  Model. 


To  meet  the  Joint  Staff’s  reporting  requirements,  the  aircraft’s  contribution 
to  mission  success  is  broken  down  one  further  level  to  a  series  of  three  components: 
the  avionics  system,  the  weapon  release  system,  and  the  aircraft  vehicle.  The  mis¬ 
sion  contribution  of  a  component  is  declared  successful  if  that  component  functions 
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properly  until  the  weapon  is  deployed.  At  that  point,  the  contribution  of  the  air¬ 
craft  to  mission  success  ends,  and  whatever  happens  to  it  after  missile  launch  is  not 
considered  in  the  model. 

In  contrast  to  the  aircraft’s  contribution,  the  weapon’s  contribution  is  broken 
down  into  two  mission  phases:  captive  carry  and  free  flight.  (The  warhead  contri¬ 
bution  is  considered  separately.)  The  weapon  must  function  properly  during  both 
phases  for  the  mission  to  be  successful.  The  captive-carry  phase  starts  when  aircraft 
power  is  applied  to  the  aircraft  avionics  during  the  operations  preflight  inspection  and 
lasts  until  the  aircraft  launches  the  missile.  The  free-flight  phase  begins  immediately 
upon  missile  launch  and  lasts  until  the  planned  warhead  detonation  time.  During 
the  free-flight  phase,  the  weapon  contribution  is  further  decomposed  into  a  series  of 
components,  each  of  which  must  function  for  the  entire  phase  to  guarantee  mission 
success.  These  components  are  the  engine,  the  fuel  system,  the  navigation/guidance 
system,  the  air  vehicle,  and  the  fuze  mechanism. 

When  the  contributions  of  all  the  delivery  systems  are  decomposed  to  the  lowest 
level,  the  system  devolves  to  a  series  structure  with  ten  independent  components. 
This  expanded  WSR  model  is  shown  in  Figure  2. 

Data  Source.  Since  it  is  obviously  not  feasible  to  collect  failure  data  under 
combat  conditions,  such  data  must  come  from  other  sources.  ACC’s  Follow-On 
Operational  Test  and  Evaluation  Phase  II  (FOT&E  II)  programs  axe  a  primary 
source  of  system  failure  data.  Since  an  FOT&E  II  program  uses  fully  operational 
versions  of  the  system  being  tested,  ACC  conducts  FOT&E  II  in  an  attempt  to  assess 
the  performance  of  its  strategic  systems  in  the  most  realistic  environment  possible. 
Currently,  this  program  provides  the  only  failure  data  that  can  be  used  to  estimate 
WSR  (although  ACC  is  currently  investigating  using  ground  test  data  to  supplement 
data  from  FOT&E  II  missions). 


Figure  2.  Expanded  WSR  Model. 

Estimating  WSR.  ACC  currently  uses  classical  statistical  theory  to  estimate 
WSR.  A  point  estimate  is  obtained  using  the  maximum  likelihood  estimator  (MLE), 
while  Easterling’s  Modified  Maximum  Likelihood  Integer  (MMLI)  estimator  is  used 
to  calculate  a  90  percent  lower  confidence  limit  [7]  (these  techniques  are  fully  ex¬ 
plained  in  Chapter  II). 

Statement  of  the  Problem 

Funding  for  ACC’s  FOT&E  II  program  has  been  substantially  reduced  in  recent 
years.  These  reductions  have  severely  limited  ACC’s  ability  to  gather  data  for  use 
in  estimating  WSR.  This  means  that,  for  some  newly  fielded  weapon  systems,  the 
failure  database  is  likely  to  remain  very  small  (that  is,  less  than  ten  data  points)  for 
some  time. 

In  the  cases  where  failure  data  is  scarce,  the  WSR  estimates  made  using  ACC’s 
current  approach  may  produce  spurious  results.  It  is  well  known  that  classical  point 
estimators  decrease  in  accuracy  as  sample  size  decreases  [16:2].  Furthermore,  because 
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the  MMLI  estimator  has  its  roots  in  asymptotic  theory,  it  is  not  intended  to  be  used 
with  small  data  sets  [12];  it  is  often  unusable  in  such  a  sampling  situation  because  of 
mechanical  difficulties  with  its  implementation  (see  Chapter  II  for  further  discussion). 

The  need  for  strategic  readiness  demands  meaningful  estimates  of  WSR.  Since 
the  data  sets  available  through  the  FOT&E  II  program  are  small,  such  estimates 
cannot  be  made  using  the  MLE/MMLI  approach.  Therefore,  the  problem  consists 
of  finding  one  or  more  acceptable  alternate  estimation  techniques  that  can  be  applied 
to  small  samples. 

Assumptions  and  Limitations 

The  following  assumptions  and  limitations  apply  to  this  research  effort: 

1.  The  context  of  this  inquiry  is  limited  to  weapon  systems  composed  of  a  combi¬ 
nation  of  bomber  aircraft  and  air-launched  missiles  (although  the  results  of  this 
study  can  easily  be  extended  to  delivery  scenarios  involving  gravity  weapons 
or  ballistic  missiles,  since  they  can  also  be  modeled  as  series  structures). 

2.  In  this  investigation,  all  the  system  components  are  assumed  to  have  bino¬ 
mial  failure  distributions.  Although  time  dependence  is  an  important  part  of 
the  WSR  model,  the  effect  of  estimating  WSR  using  a  mixture  of  pass-fail 
and  time-dependent  failure  distributions  is  not  explored.  The  only  system  re¬ 
liability  estimators  that  have  been  developed  for  the  mixed-distribution  case 
require  the  time-dependent  failure  distributions  to  be  exponential;  this  would 
require  gross  assumptions  about  the  mission  process.  Furthermore,  these  esti¬ 
mators  are  either  derived  using  asymptotic  approximations  or  are  impractical 
to  implement. 

3.  This  effort  focuses  on  one-sided  intervals  since  they  are  the  limits  of  interest 
in  most  practical  reliability  estimation  situations.  However,  it  is  possible  to 
derive  either  one-  or  two-sided  interval  estimates  for  system  reliability. 
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4.  Only  series  systems  are  considered  in  this  research  effort,  since  the  WSR  model 

is  a  series  structure. 

General  Approach 

This  investigation  explores  the  performance  of  several  classical  and  Bayesian 
alternatives  to  ACC’s  current  combination  of  estimators.  (Appendix  A  contains 
an  overview  of  Bayesian  statistics  for  the  benefit  of  the  reader  unfamiliar  with  the 
subject.)  The  performance  of  the  MLE  for  system  reliability  is  compared  to  that  of 
the  corresponding  Bayes  estimator.  Among  interval  estimators,  the  MMLI  estimator 
is  compared  to  three  alternatives;  the  Lindstrom-Madden  estimator  and  Bayesian 
probability  interval  estimators  constructed  using  approximate  beta  and  Bayes  Monte 
Carlo  empirical  posterior  distributions. 

The  relative  performance  of  the  point  estimators  is  assessed  by  comparing  their 
mean  square  errors.  For  the  four  interval  estimators,  the  interval  coverage  probability 
and  the  average  interval  lower  bound  are  compared.  The  baseline  values  of  these 
performance  measures  are  generated  using  Monte  Carlo  simulation. 

For  the  Monte  Carlo  comparison  study,  32  representative  systems  are  con¬ 
structed.  The  number  of  components,  the  underlying  component  failure  probabil¬ 
ities,  the  component  test  sample  sizes,  and  the  strength  and  accuracy  of  the  prior 
information  are  varied  across  the  systems.  Sample  sizes  are  limited  to  ten  or  less. 
Using  sets  of  simulated  failure  data,  the  reliability  of  each  system  is  estimated  re¬ 
peatedly,  and  performance  measures  generated  for  each  estimator.  The  sensitivity 
of  the  results  of  the  study  to  variations  in  system  characteristics  is  analyzed. 

Based  on  the  analysis  of  estimator  performance,  the  best  point  and  interval 
estimators  are  recommended  as  replacements  for  the  current  WSR  estimators.  These 
best  estimators  are  selected  based  on  their  performance  across  systems. 
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Sequence  of  Presentation 


Chapters  II  and  III  present  a  survey  of  system  reliability  estimation  methods. 
In  Chapter  II,  classical  point  and  interval  estimators  are  considered,  while  Chapter 
III  describes  Bayesian  estimators.  Each  major  system  reliability  estimation  approach 
is  described,  ajid  its  strengths  and  weaknesses  are  discussed  with  respect  to  small- 
sample  analysis. 

Chapter  IV  outlines  the  methodology  used  in  the  Monte  Carlo  comparison 
study.  The  simulation  design  is  constructed  and  the  methodology  of  the  compari¬ 
son  is  described  in  detail.  Technical  issues  associated  with  the  simulation  are  also 
discussed. 

Chapter  V  presents  the  findings  of  the  Monte  Carlo  study.  The  simulation 
results  are  reported  and  analyzed,  conclusions  are  drawn  about  method  performance, 
and  the  preferred  estimation  methods  are  identified.  Recommendations  for  further 
research  are  also  included. 


IJ.  Classical  System  Reliability  Estimation 


Overview 

This  chapter  surveys  the  key  results  in  classical  system  reliability  estimation 
theory.  In  keeping  with  the  scope  of  the  investigation,  the  treatment  is  limited  to 
the  case  of  a  series  system  of  independent  components,  each  with  a  binomial  failure 
distribution. 

The  survey  begins  with  point  estimation.  The  MLE  of  system  reliability  is 
derived,  and  its  asymptotic  variance  is  developed.  In  addition,  the  advantages  and 
disadvantages  of  using  the  MLE  with  small  samples  are  explored. 

The  discussion  of  point  estimation  is  followed  by  three  sections  on  interval 
estimators.  Buehler’s  exact  lower  confidence  limit  (LCL)  is  presented  first,  followed 
by  several  analytical  and  bootstrap  approximations  to  the  exact  limit.  In  each 
section,  the  estimators  are  described,  and  their  performance  with  respect  to  small 
samples  is  considered. 

Point  Estimation 

The  MLE  for  System  Reliability.  In  general,  the  reliability  of  a  component  is 
the  probability  that  it  will  function  for  a  specified  period  of  time.  If  the  component 
is  repeatedly  tested  over  a  fixed  operating  time,  its  reliability  can  be  modeled  with 
a  binomial  failure  distribution.  In  mathematical  terms,  if  the  fth  component  in  a 
system  has  a  binomial  failure  distribution,  its  reliability  is  defined  as 

Ri  =  l-Pi  (3) 

where  p,  is  the  probability  that  the  component  will  fail.  Since  a  series  system  will 
not  function  unless  all  of  its  components  function,  it  is  easy  to  see  that  the  reliability 
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of  that  system  must  be  the  product  of  the  component  reliabilities.  That  is, 

i=l  .=1 

where  k  is  the  number  of  components  in  the  system.  Since  the  pi  are  unknown 
parameters,  they  must  be  estimated  (which,  in  turn,  requires  that  the  /?,  and  R  be 
estimated).  MLEs  are  usually  used  to  obtain  the  point  estimates,  because  they  have 
convenient  mathematical  properties  like  asymptotic  normality  and  invariance  [9:53]. 

The  invariance  property  of  MLEs  is  useful  for  developing  a  point  estimator  for 
R.  This  property  states  that  an  arbitrary  transformation  of  an  MLE  space  is  the 
MLE  of  that  tranformation  (Mood  et  al.  formally  state  and  prove  the  invariance 
theorem  [33:285]).  Suppose  that,  for  component  i,  r,  failures  are  observed  in  n, 
component  tests.  Since  it  can  be  shown  that  the  MLE  for  p,  is  p,  =  r^/ni  [32:420],  the 
invariance  property  can  be  applied  to  Equation  (3)  to  get  the  MLE  of  the  component 
reliability: 

Ri  =  \-  Pi  (5) 

A  second  application  of  the  invariance  property  to  Equations  (4)  and  (5)  yields  the 
MLE  of  system  reliability: 

fl=nA  (6) 

»=i 

The  Reliability  Estimator ’s  Asymptotic  Variance.  The  variance  of  R  is  impor¬ 
tant  to  the  development  of  certain  interval  estimation  methods.  The  delta  method, 
which  is  based  on  the  asymptotic  normality  property  of  MLEs,  can  be  used  to  derive 
an  asymptotic  approximation  for  this  variance  [9].  Let  9  =  {^i, . . .  ,1?^}  be  the  MLE 
of  6  =  {9i, . . .  ,9k}.  Also,  define  the  arbitrary  function  ^  =  g{9).  By  the  invariance 
property,  ^  =  g{9)  is  the  MLE  of  It  can  be  shown  that  the  variance  of  I  can  be 
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asymptotically  approximated  by 


Cov[0„  6j] 


(■?) 


[9:52].  Equation  (7)  can  be  applied  to  a  reliability  context  by  letting  (  =  R,  6  =  {R,} 
and  9i  =  Ri.  This  yields 


{ft.) 


jdR' 

[dR,. 


Cov[Ri, 


(8) 


{fl.} 


Since  the  system’s  components  are  assumed  to  be  independent,  Equation  (8)  becomes 


'R 


k 

E 

t=i 


mJ 


R. 


(9) 


{R.} 


where  crj^  is  the  variance  of  /Z,.  This  variance  can  be  approximated  in  the  limit  by 


Riil  -  R.) 


n, 


(10) 


[12],  so  Equation  (9)  reduces  to 


t=i 


r. 


[ni(nj  -  n) 


(11) 


Small-Sample  Reliability  Point  Estimation.  The  advantages  of  using  MLEs 
instead  of  other  estimators  are  treated  at  length  in  many  advanced  statistics  texts 
(such  as  [16]).  These  advantages  are  such  that  MLEs  are  the  preferred  point  esti¬ 
mators  in  a  variety  of  situations.  Like  any  other  classical  estimator,  however,  the 
accuracy  of  the  estimates  produced  by  an  MLE  decreases  as  the  size  of  the  data  set 
used  in  the  estimate  decreases  [16:2].  Further,  the  approximate  normality  of  MLEs 
may  not  hold  for  small  samples,  so  the  asymptotic  variance  of  the  MLE  is  likely  to 
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be  a  poor  approximation  of  its  true  variance  [9:64],  Therefore,  the  MLE  of  system 
reliability  should  be  used  with  care  when  data  is  scarce. 

Interval  Estimation:  Buehler’s  Exact  LCL 

Because  the  point  estimator  of  series  system  reliability  is  the  product  of  the 
estimators  of  the  component  reliabilities,  it  seems  intuitive  that  a  confidence  interval 
for  R  could  be  constructed  by  taking  the  product  of  the  component  interval  bounds. 
While  this  approach  will  indeed  produce  a  confidence  interval,  it  will  not  yield  the 
shortest  possible  interval  estimate  for  R  [6].  The  construction  of  this  shortest  interval 
is  not  a  simple  task. 

In  1957,  R.  J.  Buehler  proposed  a  method  for  obtaining  the  exact  limit  of  the 
shortest  one-sided  confidence  interval  for  the  product  of  two  binomial  parameters  [6]. 
Buehler’s  work  was  generalized  and  applied  to  system  reliability  by  Lipow  [21]  and 
Steck  [42].  The  general  form  of  Buehler’s  method  is  examined  in  detail  in  several 
sources;  the  discussion  in  Crowder  et  al.  [9:189-191]  is  followed  here. 

Let  r  =  {ri,i  =  1, . . . ,  A:}  be  a  vector  of  component  failure  observations  from 
some  system  test.  Further,  let  p  =  =  1, . . . ,  fc}  be  any  possible  combination  of 

observed  failure  levels  for  the  system  other  than  that  defined  by  r.  Also,  define  an 
arbitrary  ordering  function  a;(-)  such  that  lower  values  of  u;(-)  are  more  desirable. 
Given  these  definitions,  Buehler’s  exact  100(1  —  7)  percent  LCL  for  the  reliability  of 
a  series  system  is  the  solution  to  the  following  nonlinear  program: 

k 

minimize  Ri  =  tt, 

»=i 

subject  to 

E 

n 

0  <  TTi  <  1,  f  =  l,...,ib  (12) 


Ui 


ni  "  I  (1  - -rirvi'"—' 
«=i  V 
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where  Q  =  {p:  uj{^  <  w(f)}.  Crowder  et  al.  suggest  that  a;(-)  =  fli  pi/^i  (the  MLE 
of  the  system  failure  probability)  is  an  intuitive  ordering  function  [9:191],  This  form 
of  u;(-)  does  not  always  produce  intuitive  results.  For  example,  Reiser  and  Jaegar 
find  that  the  resulting  interval  estimate  does  not  uniformly  shorten  as  the  component 
sample  sizes  increase;  nonetheless,  the  failure  probability  estimator  is  the  ordering 
function  most  commonly  used  in  the  solution  of  the  Buehler  minimization  problem 
137). 

Lipow  and  Riley  provide  tables  of  solutions  to  Equation  (12)  when  the  com¬ 
ponent  sample  sizes  are  equal  and  k  <  3  [22],  [23],  For  A:  >  3  or  large  n,  values, 
solving  Equation  (12)  may  not  be  possible,  or  may  require  computational  resources 
beyond  those  available  to  most  analysts  [9:191],  Willits  finds  that,  for  series  systems, 
the  nonlinear  optimization  program  MINOS  can  easily  be  used  to  find  the  solution 
to  Equation  (12),  provided  the  number  of  components  is  less  than  three  and  the 
test  samples  are  not  too  large;  however,  a  general  solution  algorithm  cannot  be  im¬ 
plemented  for  systems  with  four  or  more  components  because  of  the  computational 
complexity  of  the  problem  [48], 

Interval  Estimation:  Approximation  Methods 

Because  Buehler’s  exact  LCL  is  difficult  to  obtain  in  most  practical  experimen¬ 
tal  situations,  several  anal’  tical  methods  for  approximating  it  have  been  proposed. 
Some  of  the  more  significant  of  these  approaches  are  summarized  below. 

This  section  does  not  contain  an  exhaustive  catalog  of  analytical  approxima¬ 
tions.  Many  of  the  proposed  strategies  for  approximating  Buehler’s  limit  have  been 
shown  to  be  ineffective.  Only  those  estimators  that  have  been  considered  in  a  wide 
variety  of  sources  are  included.  The  reader  interested  in  a  comprehensive  treatment 
of  analytical  interval  estimation  methods  is  referred  to  Section  10.3  of  [27]. 
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The  Lindstrom- Madden  Estimator.  D.  L.  Lindstrom  and  J.  H.  Madden  pro¬ 
pose  a  conceptually  simple  method  for  approximating  the  exact  LCL  of  R.  In  the 
Lindstrom-Madden  method,  the  system  test  is  construed  as  a  single  binomial  exper¬ 
iment  where  n*  =  mini{n,}  tests  produce  r*  =  (1  —  R)n“‘  failures.  The  100(1  —  7) 
percent  LCL  of  R  can  then  be  calculated,  as  it  is  for  a  single  binomial  component, 
by  solving  for  R^  in  the  following  expression  [24]: 

e(".  j  (1  -  «£.)  W'”'"'’ =  7  (13) 


It  is  clear  that  r*  will  not  always  be  an  integer;  in  that  case,  Equation  (13) 
cannot  be  used  to  evaluate  Rl.  This  problem  can  be  avoided  by  using  the  identity 


r 


E 


( 


n 

i 


+  /%— (l-r)Vr 

r(n  -f  1)  yo  ^  ’ 

r,r  +  1) 


(14) 


where  /(•)  is  the  incomplete  beta  function  [17:214-215].  The  approximate  100(1  —7) 
percent  LCL  can  be  found  by  solving  the  following  equation  for  Rl  [24]; 

/(i?i,;n*-r*,r'-f  1)  =7  (15) 

Equation  (15)  can  only  be  used  when  r*  ^  n*;  fortunately,  this  condition  holds  in 
nearly  all  practical  circumstances,  since  r*  =  n*  if  and  only  if  .R  =  0. 

The  Lindstrom-Madden  estimator  has  both  weaknesses  and  strengths.  Ac¬ 
cording  to  Winterbottom,  Lindstrom-Madden  performance  deteriorates  eis  the  min¬ 
imum  component  sample  size  decreases  [50].  Also,  Martz  and  Duran  find  that  the 
Lindstrom-Madden  confidence  limits  can  be  far  too  conservative  in  certain  circum¬ 
stances  [29].  However,  Sudakov  finds  that  if  r*  is  an  integer,  the  Lindstrom-Madden 
estimator  produces  Buehler’s  exact  LCL.  Sudakov  also  shows  that  if  r*  is  not  an 
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integer,  Buehler’s  limit  is  bounded  by  Ri  (the  solution  to  Equation  (15))  and 
which  is  the  solution  to 


/(/?z,;n*-(r*,[r-  +  l)  =  7 


(16) 


where  [r*  is  the  greatest  integer  less  than  r*.  Thus,  Buehler’s  exact  LCL  is  always 
bounded  by  two  Lindstrom- Madden  limits  [43].  (Soms  claims  that  Sudakov’s  proof 
is  flawed,  but  that  in  most  practical  situations  the  Lindstrom-Maddeii  solutions  still 
bound  the  exact  limit  [39].) 


Interval  Estimation  Based  on  Asymptotic  Normality.  Since  R  is  a  MLE,  a 
straightforward  way  to  obtain  a  LCL  on  system  reliability  is  to  take  advantage 
of  the  fact  that  MLEs  are  asymptotically  normal  [16:43-44].  It  is  clear  that  if  R 

A  A 

is  normal  with  mean  R  and  variance  defined  by  Equation  (11),  then  {R  —  R)I(t^ 
has  an  asymptotic  standard  normal  distribution.  This  means  that  an  approximate 
100(1  —  7)  percent  LCL  can  be  estimated  by  solving  the  following  equation  for  Ri: 


Pr 


{R-R) 


<  Rl 


=  7 


(17) 


In  the  literature,  this  approach  is  often  referred  to  as  the  maximum  likelihood  (ML) 
method  (although  occasionally  it  is  called  “linearization”  [25]). 

An  obvious  weakness  of  the  ML  method  is  its  use  of  a  symmetric  distribution 
with  infinite  support  (the  normal  distribution)  to  approximate  a  random  variable 
with  an  unsymmetric,  bounded  distribution  {Rl)-  In  other  words,  it  is  possible  in 
some  situations  to  produce  a  LCL  that  lies  outside  the  unit  interval.  The  next  five 
estimators  presented  attempt  to  improve  this  situation. 
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The  MMLI  Estimator.  Easterling’s  MMLI  estimator  is  a  modification  of  the 
ML  approach  described  in  the  previous  section.  As  with  the  Lindstrom- Madden 
estimator,  Easterling  conceives  of  the  system  as  undergoing  a  single  binomial  exper¬ 
iment  rather  than  a  set  of  independent  component-level  experiments.  The  difference 
between  the  MMLI  and  the  Lindstrom-Madden  estimators  is  that  Easterling  de¬ 
rives  r*  and  n*  by  equating  the  variance  of  R  under  the  single  experiment  with  the 
asymptotic  variance  of  the  estimator  of  a  binomial  probability. 

If  R  is  construed  as  a  parameter  of  a  Binomial(n*,  R)  distribution,  and  R  is 
the  MLE  of  Equation  (10)  can  be  applied  to  get  the  variance  of  R: 

.2  m  -  R) 


Setting  the  variances  in  Equations  (11)  and  (18)  equal  yields 


«=I 


Vi 


ni(ni  -  r<)J 


R(1  -  R) 


n* 


which,  when  solved  for  n*,  produces  the  expression 


n  = 


(1--R) 


RT.U 


(19) 


(20) 


It  follows  immediately  that  r*  =  (1  —  R)n*. 

The  values  of  r*  and  n*  are  often  not  integers.  Beised  on  empirical  evidence, 
Easterling  concludes  that  rounding  both  r*  and  n*  up  to  the  next  integer  improves 
the  approximation  of  Buehler’s  limit.  Applying  the  identity  in  Equation  (14),  the 
100(1  —  7)  percent  LCL  can  be  found  by  solving 


/(/?i,-(K-rr*),(fr*  +  l))  =  7 


(21) 
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for  Rli  where  [r*  and  fn*  are  the  smallest  integers  greater  than  r*  and  n*,  respec¬ 
tively  [12].  It  should  be  noted  that  when  R  =  \,  n*  =  0  and  r*  is  indeterminate, 
so  the  MMLI  estimator  cannot  be  used.  With  respect  to  small-sample  analysis,  this 
is  a  major  weakness,  since  R  often  equals  one  when  rii  values  are  small  and  the  Ri 
values  are  reasonably  high. 

Asymptotic  Expansion.  Winterbottom  proposes  another  refinement  to  the  ML 
method  based  o  Cornish-Fisher  asymptotic  expansions.  This  approach  improves  the 
approximation  of  the  LCL  by  using  the  expansions  to  correct  for  bias  and  skewness. 
The  result  is  an  approximate  LCL  that,  for  certain  large-sample  situations,  performs 
very  well.  Winterbottom’s  approximate  100(1  —  7)  percent  LCL  is 


Ri  =  RK 


(22) 


where 


K  =  1 
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(^2  (l-n?)(n»)»  0  4 
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cr^  =  or'^/R'^,  n*  —  min,{n,}  and  Zi-.y  is  the  100(1  —  7)th  percentile  of  the  standard 
normal  distribution.  Winterbottom  notes  that  the  approximation  requires  >  0  Vt. 
The  derivation  of  Equation  (22)  is  omitted  here;  for  details,  see  [49]  and  [50]. 


The  Likelihood  Ratio  (LR)  Lower  Bound.  Madansky  proposes  an  asymptotic 
approximation  method  based  on  the  LR  test  statistic.  For  a  series  system,  A  =  R/ R 
is  the  LR  test  statistic  for  system  reliability.  Wilks  shows  that,  for  general  A,  —2  In  A 
has  a  limiting  distribution  with  one  degree  of  freedom  [47].  Using  Wilks’  result, 
Madansky  derives  the  following  100(1  —  7)  percent  LCL  for  R  [25j: 

RL  =  inf{R;-21nA<x"(l)}  (23) 
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Madansky’s  result,  while  of  historical  importance,  does  not  improve  the  gross  approx¬ 
imation  of  the  ML  method  enough  to  be  practically  useful.  In  addition.  Equation 
(23)  is  very  difficult  to  solve  on  a  computer  [49], [50]. 


The  Approximately  Optimal  (AO)  Nonrandomized  LCL.  Based  on  asymptotic 
theory,  Mann  develops  an  approximate  LCL  that  compensates  for  the  case  where 
r,-  =  0  for  some  i.  (For  a  series  system,  most  other  approximation  methods  ignore  the 
impax:t  of  this  case.)  Mann’s  derivation,  which  is  beyond  the  scope  of  this  discussion, 
centers  around  the  development  of  asymptotic  approximations  for  first  two  moments 
of  the  random  variable  —  In  R. 


If  fl  is  the  set  of  components  with  at  least  one  observed  failure  in  a  sample 
(i.e.,  fl  =  {t  ;  1  <  i  <  fc,r,  ^  0}),  then  the  approximations  of  the  first  two  moments 
of  —  In  iZ  are  given  by 


0.5[1 -f  (l/i/)]  1-^ 

N  1  -  0.5(1  -  R) 


(24) 


and 


/^2 


0.5[1  +  il/u)]pi 
N 


(25) 


where 

N  =  [inf{ni}]  [l  -  0.5(1  -  Rf]  [l  -  0.5(1  -  .ff)]  (26) 


and 


(27) 


Mann  applies  Equations  (24)  and  (25)  to  get  the  following  expression  for  the  AO 
nonrandomized  100(1  —  7)  percent  LCL  [26], [27]: 


Rl  =  exp 


(28) 
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Poisson  Approximation.  When  the  component  reliabilities  are  known  to  be 
high,  Buehler’s  exact  LCL  can  be  approximated  by  assuming  each  has  a  Poisson 
distribution.  If  the  sample  sizes  are  equal  and  large,  it  can  be  shown  that 

=  (29) 

where  i/  =  2[(X^fLi  r.)  +  1]  and  n*  =  min,{n,};  if  the  sample  sizes  differ,  Equation 
(29)  produces  a  conservative  LCL  [50].  The  Poisson  assumption  breaJcs  down  when 
samples  are  small;  for  a  detailed  treatment  of  Poisson-related  techniques,  see  Section 
10.3  of  [27]. 

(/sing  Approximation  Methods  with  Small  Samples.  All  the  interval  estimators 
presented  in  this  section,  except  the  Lindstrom-Madden  estimator,  are  based  on 
asymptotic  theory.  Asymptotically-derived  estimators  are  meant  to  be  used  with 
sufficiently  large  data  sets;  when  sample  sizes  are  small,  the  approximation  is  often 
inaccurate  [9:64].  In  some  cases,  as  with  the  Poisson  approximation,  the  inaccuracy 
is  severe;  in  other  instances,  as  with  ML  and  LR,  the  approximations  are  not  the 
best  available  even  when  they  are  used  with  a  large  sample  [50]. 

Some  approximate  estimators  suffer  because  they  fail  to  account  for  the  case 
where  no  failures  are  observed  for  one  or  more  components  [26].  When  r,  =  0,  ML, 
MMLI  and  LR  all  ignore  the  contribution  of  component  i  to  Ri,\  when  r,  =  0  for 
every  component,  these  estimators  cannot  be  used  at  all.  Also,  asymptotic  expansion 
is  useless  if  r,  =  0  for  one  or  more  components.  Mann’s  AO  LCL  is  an  attempt  to 
rectify  this  situation,  but  its  derivation  is  still  based  on  asymptotic  theory. 

The  only  approximate  estimator  that  can  be  legitimately  used  for  small  data 
sets  is  the  Lindstrom-Madden  estimator.  While  this  method  has  its  drawbacks,  it  is 
nonetheless  a  viable  alternative  to  the  MMLI  method  for  small  test  samples. 
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Interval  Estimation:  Bootstrap  Methods 


According  to  Crowder  et  al.,  the  use  of  Monte  Carlo  simulation  to  approximate 
system  reliability  confidence  limits  dates  to  the  1963  work  of  Rosenblatt.  Efron 
refined  and  extended  the  method,  and  is  generally  credited  with  originating  the 
term  “bootstrap”  [9:192].  In  bootstrap  interval  estimation,  Monte  Carlo  random 
variates  are  generated  and  used  to  construct  an  empirical  cumulative  distribution 
function  (CDF).  The  appropriate  quantile  of  this  CDF  is  used  to  approximate  the 
exact  confidence  limit.  The  steps  in  the  general  procedure,  as  implemented  by  Levy 
and  Moore,  are  as  follows: 

1.  Determine  the  underlying  component  failure  distributions  from  the  observed 
data  set  using  the  methods  of  ML. 

2.  Generate  a  set  of  simulation  component  reliabilities  using  the  component  dis¬ 
tributions  identified  in  Step  1. 

3.  Calculate  the  simulated  system  reliability  as  a  function  of  the  component  reli¬ 
abilities. 

4.  Repeat  Steps  2  and  3  iV  —  1  times,  where  N  is  some  large  number  (say,  10000). 

5.  Rank  order  the  simulated  system  reliabilities,  with  the  lowest  magnitude  first. 

6.  Use  the  (iV7)th  simulated  system  reliability  as  the  approximate  100(1  —  7) 
percent  LCL. 

While  Levy  and  Moore  are  concerned  primarily  with  continuous  component  fail¬ 
ure  distributions,  the  above  procedure  is  applicable  to  systems  with  binomially- 
distributed  components  [20]. 

Asymptotic  Normal  Monte  Carlo  (ANMC).  Rice  and  Moore  proposed  a  vari¬ 
ant  of  the  bootstrap  method  for  the  all-binomial  case.  This  refinement  recognizes 
that  since  the  parameters  of  the  underlying  failure  distributions  are  MLEs,  their 
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values  can  vary  under  repeated  sampling.  The  asymptotic  normality  of  these  MLEs 
is  used  to  account  for  the  estimation  uncertainty.  (The  refined  method  also  compen¬ 
sates  for  cases  when  r,  =  0.)  The  steps  in  ANMC  are  as  follows: 

1.  From  the  observed  test  data,  determine  for  i  =  1, . . . ,  fc.  If  =  0  for  some 
i,  calculate  r*  by  solving 

(r*)^  +  (2n.  -  3)r*  -j-  (n,-  -  1)  ln[X(i_-,;(  ’  '1  =  0-  (30) 


2.  Calculate 


R.  =  { 


1  -  in/ni)  if  r,  >  0 
[  1  -  {r^/ni)  if  r,  =  0 


and 


.2 

O:  —  - . 

m 


(31) 


(32) 


3.  Set  N  large  (say,  10000).  For  j  =  1, . . . ,  iV,  generate  a  set  of  k  standard  normal 
deviates  Zi,j  (t  =  l,...,k). 

4.  Calculate  -f  Ri,  V(j,  j). 

5.  Calculate  Rj  =  n?=i  V(i,  j). 

6.  Rank  order  the  values  of  Rj  to  get  the  order  statistics  R^jy 

7.  Approximate  the  lower  limit  as  Ri  =  R(/V7)- 


Rice  and  Moore  observe  a  slight  upward  bias  in  the  confidence  intervals  ANMC 
produces  [38]. 


The  Modified  ANMC  Method.  Chao  and  Huwang  propose  a  modification  to 
ANMC  to  correct  its  inherent  bias.  In  modified  ANMC,  Bayes  estimates  are  used 
for  the  component  reliabilities.  (Although  this  makes  the  modified  method  pseudo- 
Bayesian,  it  is  considered  here  instead  of  in  Chapter  III  because  it  is  intended  to  be 
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used  in  a  classical  context.)  The  steps  for  the  modified  method  are  the  same  as  those 
for  the  original  method,  except  that  the  following  terms  are  defined  differently: 


R,  =  l- 


ri  + a 
Ui  +  Q  +  /3 


=  1  - 


Pi  +  a 
Mi  +  a  +  /? 


where  pi  ~  Bi(ni,^),  and  a  and  are  fixed  prior  parameters.  Based  on  empirical 
evidence,  Chao  and  Huwang  find  that  the  modified  ANMC  method  produces  the 
best  interval  coverage  if  0.1  <  a  <  0.2  and  =  0  [8]. 


Using  the  Bootstrap  with  Small  Samples.  Martz  and  Duran  find  that  the  per¬ 
formance  of  the  standard  bootstrap  varies  widely,  depending  on  the  situation  in 
which  it  is  used;  in  many  cases,  it  produces  quite  optimistic  intervals  [29].  Also,  the 
standard  bootstrap  does  not  accomodate  the  case  where  one  or  more  components  has 
no  failures  observed;  if  r,  =  0  for  some  i,  the  component  is  treated  as  if  it  had  per¬ 
fect  reliability.  Rice  and  Moore  rectified  this  situation  with  the  ANMC  method  [38]. 
Unfortunately,  the  device  they  use  to  incorporate  the  variability  of  the  reliability 
estimates  is  based  on  the  asymptotic  normality  of  the  estimators,  and  their  method 
produces  slightly  optimistic  LCLs.  Chao  and  Huwang  avoid  the  use  of  asymptotic 
theory  when  correcting  the  bias;  for  component  sample  sizes  of  ten,  the  method 
appears  to  yield  reasonable  results  [8].  However,  no  computational  experience  has 
been  reported  for  smaller  samples  or  for  unequal  sample  sizes. 
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III.  Bayesian  System  Reliability  Estimation 


Overview 

This  chapter  summarizes  Bayesian  system  reliability  estimation  theory.  As 
with  Chapter  II,  the  methods  considered  are  limited  to  the  case  of  a  series  system 
of  independent  binomial  components.  After  introductory  remarks  concerning  the 
selection  of  the  component  prior  densities  and  the  nature  of  reliability  estimation, 
methods  for  determining  the  posterior  density  of  R  are  surveyed  and  evaluated.  The 
chapter  concludes  with  a  brief  justification  for  investigating  Bayesian  alternatives  to 
classical  reliability  techniques. 


Selecting  the  Prior  Density  of  R 

When  a  component  has  an  underlying  binomial  failure  density,  one  possible 
prior  probability  density  function  (pdf)  is  the  Beta(a,/?)  density 


9{v)  = 


r(a  +  ^) 

ria)m 


„ar-l 


(1  -p) 


(35) 


(other  forms,  such  as  the  uniform  density,  are  used  as  well).  Since  a  beta  prior  is 
a  natural  conjugate  prior  for  a  binomial  proportion,  it  has  desirable  mathematical 
properties  (see  Appendix  A).  Further,  Weiler  showed  that  assumption  of  a  beta 
prior  when  the  true  form  is  other  than  beta  has  negligible  effect  on  the  posterior 
of  R  [46].  Because  of  these  advantanges,  the  beta  family  is  most  often  used  in  a 
binomial  sampling  situation. 

When  the  Beta(a,/9)  density  is  used  as  the  prior  pdf  of  a  component,  its  hy¬ 
perparameters  are  usually  chosen  as  a  =  s®  -1- 1  and  13  =  n°  —  s°  +  1.  The  values  of 
s°  and  n°,  which  represent  the  number  of  successes  and  trials  in  a  prior  pseudo-test 
of  component  i,  attempt  to  adequately  translate  the  subjective  prior  information 
about  the  reliability  of  component  i  into  numerical  values  [30:236].  The  values  of 
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3°  and  n°  can  be  determined  somewhat  objectively  (using,  for  example,  test  data 
from  a  similar  system  or  component),  or  the  estimation  process  may  be  completely 
subjective. 

If  the  subjective  approach  must  be  used,  the  prior  hyperparameters  should  not 
be  assigned  values  directly.  Mosleh  and  Apostolakis  recommend  that  the  source  of 
the  prior  assessment  be  the  perceived  quantiles  of  the  prior  distribution.  This  is 
because  the  resulting  prior  variance  tends  to  be  too  small  otherwise  [34].  Building 
on  research  by  Weiler  [46],  Waterman  et  al.  (reproduced  in  [30:236-237])  propose 
an  assessment  method  for  beta  hyperparameters  which  uses  the  mean  and  one  other 
quantile.  The  procedure  requires  that  the  assessor  assign  perceived  values  of  the 
mean  and  the  fifth  or  95th  percentile  of  the  prior  density.  Using  tables  provided  by 
Waterman  et  al.,  these  assignments  can  be  translated  into  equivalent  values  of  s° 
and  n°. 

Using  Waterman’s  procedure,  Duran  and  Booker  examine  the  effect  on  system 
reliability  estimation  of  faulty  subjective  assessments  of  beta  prior  hyperparameters. 
They  find,  not  surprisingly,  that  as  the  accuracy  of  the  prior  assessment  worsens, 
the  posterior  variance  of  R  increases  [11].  Since  this  adversely  affects  the  quality 
of  the  system  reliability  estimates,  extreme  care  should  be  taken  when  subjectively 
assessing  the  hyperparameters  of  the  component  prior  densities. 

Bayes  Estimators  and  Probability  Intervals 

When  the  posterior  pdf  of  R  is  determined,  point  and  interval  estimates  for 
R  can  easily  be  derived.  As  was  shown  in  Appendix  A,  the  Bayes  estimator  is  the 
mean  of  the  posterior  pdf  if  a  squared-error  loss  function  is  assumed.  Similarly,  the 
100(1  —  7)  percent  Bayesian  lower  probability  limit  is  the  (1007)th  percentile  of 
the  posterior  CDF,  which  is  found  by  solving  Equation  (49)  (see  Appendix  A).  Since 
these  calculations  are  usually  straightforward  when  the  form  of  the  posterior  density 
is  known,  the  primary  task  in  Bayesian  system  reliability  estimation  is  to  derive  the 


23 


posterior  for  R.  Three  methods  for  generating  this  posterior  can  be  applied  to  the 
all-binomial  ca^e.  These  methods  are  described  in  the  next  section. 


Methods  for  Deriving  the  Posterior  Density  of  R 

The  Exact  Posterior.  Springer  and  Thompson  propose  an  analytical  method 
for  deriving  the  posterior  distribution  of  series  system  reliability  in  the  case  where 
all  components  have  binomial  failure  distributions.  Since  the  posterior  derived  using 
this  method  is  exact,  the  Bayesian  point  and  interval  estimates  it  yields  are  also 
exact . 

Melhn  transform  theory  can  be  readily  applied  to  the  problem  of  obtaining  the 
posterior  density  for  R  (see  Appendix  B  for  a  discussion  of  the  Mellin  transform).  Let 
M{f{x)\u]  be  the  Mellin  transform  of  /(x)  with  respect  to  the  complex  parameter 
u.  Springer  and  Thompson  show  that 

M{g{R\s)-,u}  =  M{gi{Ri\siy,u}  (36) 

«=i 


where  s  =  {sj, . . . ,  s/t}.  If  /Z,  has  a  Beta(s°-j-l,  n"  — s°-l-l)  prior,  the  Mellin  transform 
of  the  posterior  pdf  of  R  is 


M{g{R\s)-,u} 


n 


(n-  +  1)! 


_ 1 _ 

n“r5(u  +  s  +  j)'> 


(37) 


where  s*  =  s°  -f-  s,-,  n*  =  n°  -|-  n,,  s  =  max.  s*,  n  =  max,  n*,  Re(u)  >  —  s,  and  Cj 
is  the  number  of  times  (u  -|-  s  -1- j)  occurs  over  the  component  Mellin  transforms. 
The  posterior  pdf  of  R,  which  is  not  stated  here,  can  be  obtained  by  evaluating  the 
inversion  integral  of  the  Mellin  transform  in  Equation  (37). 

The  mean  of  the  posterior  pdf  is 


£[/Z|s]  =  M{giR\s;u  = 


(38) 
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This  expression  is  also  the  Bayes  estimator  of  R  for  a  squared-error  loss  function 
[41]. 

The  closed  form  of  the  posterior  CDF  of  R  is  obtained  by  integrating  the  pdf. 
Springer  and  Thompson  show  that  this  CDF  is 


j=i 


i=l 


(39) 


where 


and 


1 

=  a  _  1  p  ^ 


t+m 


O' - 1)!  ,;^o 


(40) 


—  \ 


(Cj-i-my.  ^  +  jpM{g{R\s)-,  u}]|u=-,+j 

0  Cj  =  0 


(41) 


Springer  and  Thompson  develop  computer  code  to  evaluate  quantiles  of  this  CDF. 
However,  their  experience  (and  that  of  others,  notably  Martz  et  al.  [31])  shows  that 
such  implementations  are  numerically  unstable,  especially  for  systems  with  high 
reliability. 

In  addition  to  the  application  of  Mellin  transforms,  other  methods  exist  for 
finding  the  exact  posterior  CDF  of  R.  While  these  approaches  are  of  theoretical 
interest,  they  are  not  without  their  difficulties.  Tang  and  Gupta  provide  a  set  of 
recurrence  relations  that  can  be  used  to  find  the  exact  pdf  of  a  product  of  indepen¬ 
dent  beta  random  variables  [44];  the  problem  with  this  approach  is  that  it  requires  a 
cumbersome  term-by-term  integration  of  the  resulting  pdf  to  find  the  CDF.  As  an¬ 
other  alternative,  Wolfe  uses  a  probabilistic  argument  to  develop  a  multiple  integral 
expression  for  the  exact  posterior  CDF  [51];  however,  when  k  is  large,  the  dimen¬ 
sionality  of  the  integral  makes  it  difficult  to  solve.  A  detailed  discussion  of  these 
methods  is  omitted.  The  interested  reader  is  directed  to  the  appropriate  references. 
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The  Approximate  Beta  Posterior.  Computationally,  it  is  very  difficult  to  find 
the  exact  posterior  of  R.  Springer  and  others  suggest  overcoming  this  difficulty  by 
approximating  the  exact  posterior  density  of  R  with  a  beta  density  sharing  the  same 
first  two  moments  [40:268-269], [45], [50].  According  to  Springer,  such  approximations 
can  often  be  quite  close  to  the  exact  form  [40:268].  Thompson  and  Haynes  report 
that,  at  least  for  moderate  sample  sizes  (n^  >  20),  this  beta  approximation  yields 
excellent  results  [45]. 

Using  Springer’s  approach,  the  approximate  beta  posterior  can  be  constructed 
easily.  The  first  two  moments  of  the  exact  posterior  pdf  of  R  can  be  found  by  setting 
w  =  2, 3  in  the  density’s  Mellin  transform.  The  resulting  moments  are 


=nf^) 

MU +  2/ 


(42) 


and 


f  A  (^‘  +  1)W  +  2)  ] 
lM«  +  2)(n:  +  3)J 


(43) 


-|-2)(nr-f  3) 

The  Beta(a,  0)  density  with  these  same  first  and  second  moments  has  the  parameters 


Q  = 


-  It)- (r^  It 


(44) 


and 


[40:269]. 


(9  = 


p{l  -  -  Or^l  -  fl) 


(45) 


The  Bayes  Monte  Carlo  Empirical  Posterior.  Martz  and  Duran  propose  a 
Monte  Carlo  method  for  approximating  a  100(1  —  7)  lower  probability  limit.  This 
approach,  which  is  a  variant  of  the  bootstrap,  has  as  its  goal  the  generation  of  an 
empirical  posterior  CDF  for  R.  While  Martz  and  Duran  use  noninformative  priors 
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in  their  implementation,  their  method  can  be  easily  adapted  to  handle  informative 
priors. 

The  steps  in  Bayes  Monte  Carlo  are  as  follows: 

1.  Set  N  large  (say,  10000).  For  j  =  generate  the  random  variate  Rij 

from  a  Beta(n,  — r,4-ai,rj+/9i)  distribution,  where  a,  and  are  the  parameters 
of  the  ith  prior  density. 

2.  Calculate  Rj  =  n?=i  Ri,ji  Vj. 

3.  Rank  order  the  values  of  Rj  to  obtain  the  order  statistics  R^jy 

4.  Approximate  the  100(1  —  7)  lower  probability  limit  as  Rl  —  R^Ny)  [29]. 

Martz  and  Duran  find  that  this  method  generally  outperforms  the  classical  bootstrap 
[29]. 

Bayes  Monte  Carlo  is  only  intended  to  produce  approximate  probability  limits. 
If  a  Bayes  estimate  of  R  is  desired.  Equation  (38)  (the  first  moment  of  Springer  and 
Thompson’s  exact  posterior  distribution)  should  be  used. 

The  Bayesian  Estimation  of  WSR 

Martz  and  Waller  argue  convincingly  for  the  use  of  Bayesian  inference  in  relia¬ 
bility  analysis.  They  claim  that  the  Bayesian  approach  has  the  following  advajitages 
over  classical  reliability  techniques: 

1.  If  the  prior  information  is  valid,  the  Bayesian  inferences  are  more  accurate. 

2.  Bayesian  analyses  usually  require  smaller  samples. 

3.  In  a  Bayesian  test,  unacceptable  results  are  a  consequence  of  bad  assumptions 
(i.e.,  inaccurate  prior  information)  rather  than  questionable  methodology. 

4.  Bayesian  probability  intervals  for  R  are  easier  to  obtain  and  understand  than 
classical  confidence  intervals  [30:170-174]. 
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On  the  face  of  it,  Martz  and  Waller’s  arguments  are  persuasive,  but  the  Bayesian 
approach  is  far  from  universally  accepted.  In  fact,  controversy  has  raged  between 
Bayesian  and  classical  statisticians  for  many  years  as  to  the  validity  of  the  list  of 
“advantages”  given  above  (see  Appendix  A). 

The  discussion  in  Chapter  II  showed,  however,  that  (1)  ML  estimation  perfor¬ 
mance  suffers  when  sample  sizes  decrease,  and  (2)  the  Lindstrom-Madden  estimator 
provides  the  only  classical  interval  estimator  that  can  be  legitimately  applied  when 
component  test  samples  are  small.  In  the  face  of  these  difficulties,  it  seems  reasonable 
to  accept  the  Bayesian  paradigm  as  a  valid  alternative  to  frequentist  theory.  This 
increases  the  list  of  potential  replacement  estimators,  and  avoids  having  to  accept 
the  MLE/ Lindstrom-Madden  combination  by  default. 

The  problem  is  that  there  is  no  small-sample  comparison  between  Bayesian 
and  classical  system  reliability  estimation  techniques  published  in  the  reliability  lit¬ 
erature.  The  remainder  of  this  thesis  is  devoted  to  conducting  such  a  study. 
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IV.  The  Monte  Carlo  Comparison  Study 


Overview 

This  chapter  describes  the  Monte  Carlo  study  undertaken  in  this  research  effort. 
First,  the  estimation  methods  that  were  compared  are  identified,  as  are  the  measures 
of  merit  upon  which  the  comparison  was  based.  The  series  systems  used  in  the  study 
axe  described  next,  followed  by  a  discussion  of  technical  considerations  related  to 
the  Monte  Carlo  implementation.  The  chapter  concludes  with  a  description  of  the 
methodology  used  in  the  simulation. 

Candidate  Estimation  Methods 

The  purpose  of  this  study  was  to  find  some  set  (or  sets)  of  point  and  upper-tail 
interval  estimators  that  performs  better  for  small  test  samples  than  those  presently 
used  by  ACC.  Accordingly,  the  performance  of  several  classical  and  Bayesian  estima¬ 
tors  were  considered  as  possible  replacements  for  both  the  MLE  (R)  and  the  MMLI 
interval  estimator. 

The  Bayes  estimator  in  Equation  (38)  was  examined  as  an  alternative  to  the 
MLE.  This  estimator  was  chosen  for  comparison  because  it  is  the  first  momei  .  of  all 
three  posterior  forms  considered  in  Chapter  III. 

Two  forms  of  the  Bayesian  interval  estimator  in  Equation  (49)  (see  Appendix 
A)  were  compared  to  the  MMLI  estimator.  The  first  form  was  constructed  from  the 
approximate  beta  posterior  CDF;  the  second  used  the  BMC  empirical  posterior  CDF 
(the  exact  posterior  was  not  used  because  of  the  computational  difficulties  associated 
with  its  implementation).  Since  the  Lindstrom-Madden  estimator  is  the  only  classical 
interval  estimator  suitable  for  use  with  small  test  samples,  it  was  also  compared  to 
the  MMLI  estimator.  All  interval  estimators  considered  were  constructed  at  two 
levels  of  significance  (0.10  and  0.05)  to  aid  sensitivity  analysis. 
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Measures  of  Merit 

The  measure  of  merit  used  to  compare  the  point  estimators  was  mean  square 
error  (MSE)  (see  Appendix  A).  Since  the  Bayes  estimator  under  consideration  min¬ 
imizes  the  squared-error  loss  function,  MSE  was  the  logical  measure  to  use.  The 
point  estimator  that  produced  the  smallest  estimated  MSE  for  a  given  system  was 
considered  the  superior  performer  for  that  system. 

For  the  interval  estimators,  the  measures  of  merit  were  the  coverage  propor¬ 
tion  (the  percentage  of  generated  intervals  containing  the  true  system  reliability) 
and  the  average  interval  lower  bound.  These  measures  are  often  used  to  compare 
interval  estimators  (see,  for  example,  Martz  and  Duran  [29]).  An  interval  estimator 
that  produced  a  coverage  proportion  less  than  (1  —  7)  was  considered  optimistic. 
For  a  given  system,  the  interval  estimator  that  produced  the  largest  average  inter¬ 
val  bound  at  both  significance  levels  without  being  optimistic  was  considered  the 
superior  estimator  for  that  system. 

Hypothetical  System  Structures 

The  basic  concepts  of  experimental  design  were  used  to  construct  a  comprehen¬ 
sive  set  of  hypothetical  series  systems  for  study.  Table  1  shows  the  design  strategy 
that  was  used  to  develop  the  systems.  This  strategy  was  based  on  a  completely 
randomized  2®  factorial  design.  Although  a  more  complex  approach  could  have  been 
taken,  only  two  factor  levels  were  considered  so  that  the  number  of  systems  remained 
manageable.  Table  2  presents  the  details  of  the  randomized  design. 

The  hypothetical  systems  were  developed  using  the  components  described  in 
Table  3.  These  components  were  constructed  so  that  their  true  reliabilities  were 
between  0.90  and  0.999,  a  reasonable  range  for  many  real-world  applications.  The 
strength  and  accuracy  of  the  component  prior  information  were  determined  by  the 
values  chosen  for  s°  and  n°.  For  the  purposes  of  this  study,  it  was  assumed  that 
setting  n°  ss  25  provided  sufficiently  strong  prior  information;  for  a  weak  prior,  n° 
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was  fixed  at  around  10.  The  accuracy  of  the  prior  information  was  controlled  by 
varying  the  magnitude  of  the  difference  between  is°/n°)  and  the  true  component 
reliability  Hi  (that  is,  |(s°/n°)  —  To  provide  a  component  with  accurate  prior 
information,  s°  and  n°  were  selected  so  that  |(sf/nf)  —  <  0.2;  an  inaccurate 

prior  was  constructed  by  choosing  the  hyperparameters  |(s°/n°)  —  i2,|  >  0.2.  Table  4 
shows  the  component-by-component  structure  of  each  series  system. 


Table  1.  Design  Strategy. 


Factor 

Level  1 

Level  2 

k 

5 

10 

Sample  sizes 

n,-  =  4  V  i 

m  =  10  V  i 

Ri 

mmmmmrn 

unequal 

Component  Prior 
Info  Strength 

uniformly  weak 

uniformly  strong 

Component  Prior 
Info  Accuracy 

uniformly  inaccurate 

uniformly  accurate 
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Table  2.  Design  Strategy  Implementation. 


System 

k 

Hi 

ft. 

Component  Prior  Information 

Strength 

Accuracy 

1 

5 

a 

weak 

inaccurate 

2 

5 

4 

equal 

weak 

accurate 

3 

5 

4 

equal 

strong 

inaccurate 

4 

5 

4 

equal 

strong 

accurate 

5 

5 

4 

unequal 

weak 

inaccurate 

6 

5 

4 

unequal 

weak 

accurate 

7 

5 

4 

unequal 

strong 

inaccurate 

8 

5 

4 

unequal 

strong 

accurate 

9 

5 

10 

equal 

weak 

inaccurate 

10 

5 

10 

equal 

weak 

accurate 

11 

5 

10 

equal 

strong 

inaccurate 

12 

5 

10 

1 

equal 

strong 

accurate 

13 

5 

10 

unequal 

weak 

inaccurate 

14 

5 

10 

unequal 

weak 

accurate 

15 

5 

10 

unequal 

strong 

inaccurate 

16 

5 

10 

unequal 

strong 

accurate 

17 

10 

4 

equal 

weak 

inaccurate 

18 

10 

4 

equal 

weak 

accurate 

19 

10 

4 

equal 

strong 

inaccurate 

20 

10 

4 

equal 

strong 

accurate 

21 

10 

11 

weak 

inaccurate 

22 

10 

4 

weak 

accurate 

mmi 

10 

Kl 

unequal 

strong 

inaccurate 

24 

10 

4 

unequal 

strong 

accurate 

25 

10 

10 

equal 

i  weak 

inaccurate 

26 

10 

10 

equal 

weak 

accurate 

27 

10 

10 

equal 

,  strong 

inaccurate 

28 

10 

10 

equal 

i  strong 

accurate 

29 

10 

10 

unequal 

weak 

inaccurate 

10 

10 

unequal 

weak 

accurate 

31 

10 

10 

unequal 

strong 

inaccurate 

32 

10 

10 

unequal 

strong 

accurate 

32 
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ft 


0.9960 

0.9510 

0.7738 

0.5905 

0.9292 

0.7237 

0.8106 

0.7931 

0.5905 

0.7738 

0.9510 

0.9960 

0.9684 

0.8757 

0.6579 

0.8448 

0.9044 

0.9900 

0.5987 

0.3487 

0.7072 

0.5105 

0.9548 

0.5293 

0.5987 

0.3487 

0.9900 

0.9044 

0.6347 

0.5920 

0.6004 

0.8211 


Technical  Considerations 

The  Monte  Carlo  software  was  written  in  FORTRAN  77;  the  accuracy  of  each 
module’s  output  was  verified  using  a  set  of  hand-computed  results.  Most  of  the 
statistical  calculations  required  by  the  simulation  were  performed  using  the  Interna¬ 
tional  Mathematical  and  Statistical  Libraries  (IMSL)  Stat/Library.  IMSL  routines 
called  by  the  simulation  software  included  the  binomial  deviate  generator  DRNBIN 
[14:967],  the  beta  deviate  generator  DRNBET  [14:993],  and  the  inverse  beta  function 
DBETIN  [14:915]. 

IMSL  Stat/Library ’s  uniform  random  number  stream,  which  it  uses  to  generate 
other  random  deviates,  was  evaluated  for  independence  and  uniformity.  A  stream 
of  random  numbers  was  generated  on  each  computer  used  in  the  simulation;  each 
stream  was  then  evaluated  using  a  runs-up  test  and  a  chi-square  goodness-of-fit  test. 
Each  test  was  conducted  using  the  guidelines  in  Section  7.4  of  Law  and  Kelton  [19]. 
The  tests  showed  that  the  uniformity  and  independence  assumptions  were  justified. 

Methodology 

The  following  series  of  steps  was  executed  for  each  of  the  32  hypothetical 
systems: 

1.  A  set  of  system  failure  data  was  constructed  by  simulating  a  binomially- 
distributed  failure  count  for  each  of  the  system’s  k  components.  Five  thousand 
of  these  data  sets  were  generated. 

2.  ML  and  Bayes  estimates  of  R  were  calculated  using  each  of  the  5000  data  sets; 
the  MSE  was  then  calculated  for  each  estimate. 

3.  MMLI  90  and  95  percent  lower  bounds  on  R  were  computed  using  each  of  the 
5000  data  sets.  Each  data  set  with  no  observed  failures  for  any  component 

A 

(that  is,  where  R  =  1.0)  was  ignored.  If  more  than  500  data  sets  were  ignored, 
the  coverage  proportion  and  the  average  lower  bound  were  not  reported;  other- 
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wise,  these  measures  were  calculated  using  only  those  lower  bounds  that  were 
evaluated. 

4.  Lindstrom-Madden  90  and  95  percent  lower  bounds  on  R  were  computed  using 
each  of  the  5000  data  sets,  and  the  coverage  proportion  and  the  average  lower 
bound  were  calculated. 

5.  The  performance  measures  for  the  approximate  beta  posterior  90  and  95  per¬ 
cent  lower  proportion  limits  on  R  were  derived  using  the  methods  in  Step  4. 

6.  For  each  of  the  5000  data  sets,  a  BMC  empirical  posterior  CDF  of  R  was 
constructed  according  to  the  method  described  in  Chapter  III.  Each  CDF  was 
composed  of  10000  data  points.  The  appropriate  order  statistics  were  used  as 
the  90  and  95  percent  lower  proportion  bounds  on  R. 

Using  these  six  steps,  a  set  of  performance  measures  was  compiled  for  each  estimator 
under  consideration.  These  performance  measures  are  reported  and  analyzed  in  the 
next  chapter. 
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V.  Findings 


Overview 

This  chapter  presents  the  findings  of  the  Monte  Carlo  comparison  study.  In 
the  first  section,  the  performance  measures  of  the  point  estimators  are  reported  and 
analyzed;  the  same  is  done  in  the  second  section  for  the  interval  estimators.  After 
the  analysis  of  the  results,  conclusions  are  drawn  about  estimator  performance.  The 
chapter  ends  with  a  list  of  recommendations  for  further  research. 

Results  and  Analysis:  Point  Estimators 

Table  5  lists  the  MSE  by  system  for  the  ML  and  Bayes  estimators  of  R.  The 
ratio  of  the  losses  (ML  to  Bayes)  is  also  included;  if  the  ratio  is  greater  than  one, 
the  Bayes  loss  is  smaller  than  that  of  the  MLE  (and  therefore  the  Bayes  estimator 
performed  better).  The  MLE  outperformed  the  Bayes  estimator  for  23  of  the  32 
systems.  Appendix  C  contains  a  discussion  of  the  Monte  Carlo  sampling  statistics 
for  the  point  estimators. 

An  analysis  was  conducted  to  determine  if  any  of  the  structural  aspects  of 
the  hypothetical  systems  had  a  statistically  significant  effect  on  the  choice  of  the 
best  estimator  for  a  given  system.  Contingency  tables  were  used  because  of  the 
nominal  nature  of  the  explanatory  variables  (the  structural  factors  in  Table  1)  and 
the  response  variable  (the  estimator  judged  best).  The  analysis,  which  followed  the 
guidelines  in  Chapters  2  and  3  of  Agresti  [1],  was  conducted  at  significance  level  0.05. 
Only  main  factor  effects  and  pair-wise  factor  interactions  were  investigated;  higher 
level  interactions  were  assumed  to  be  negligible.  The  results  of  the  analysis  are  in 
Tables  6  and  7;  the  relevant  contingency  tables  «ire  in  Appendix  D. 
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Table  5.  Point  Estimator  Performance  Measures. 


System 

- 1 

MSEiR) 

— 

MSE(Bayes) 

Ratio  of  MSEs 
(ML  to  Bayco) 

1 

0.0012 

0.5747 

0.0021 

2 

0.0114 

0.0562 

0.2026 

3 

0.0414 

0.2573 

0.1608 

:  4 

0.0504 

0.0162 

3.1052 

5 

0.0158 

0.4867 

0.0324 

6 

0.0463 

0.0197 

2.3504 

7 

0.0374 

0.3060 

0.1224 

8 

0.0393 

2.3447 

9 

0.0203 

0.1506 

0.1347 

10 

0.0160 

0.0079 

2.0290 

11 

0.0046 

0.3676 

0.0125 

12 

0.0005 

0.0155 

0.0307 

13 

0.0030 

0.3695 

0.0081 

14 

0.0102 

0.0190 

0.5366 

15 

0.0197 

0.1467 

0.1345 

16 

0.0126 

0.0116 

1.0830 

17 

0.0211 

0.7224 

0.0292 

18 

0.0025 

0.2181 

0.0.113 

19 

0.0375 

0.0938 

0.3995 

20 

0.0509 

0.0149 

3.4237 

21 

0.0471 

0.4424 

0.1065 

22 

0.0499 

0.0427 

1.1699 

23 

0.0105 

0.7710 

0.0136 

24 

0.0491 

0.0326 

1.5062 

25 

0.0194 

0.2474 

0.0785 

26 

0.0138 

0.0117 

1.1712 

27 

0.0009 

0.7528 

0.0012 

28 

0.0106 

0.0191 

0.5543 

29 

0.0201 

0.3052 

0.0658 

30 

0.0196 

0.0282 

0.6964 

31 

0.0184 

0.3388 

0.0544 

32 

0.0136 

0.0262 

0.5204 
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Table  6.  Main  Factor  Effects:  Point  Estimators. 


Factor 

P-Value 

Significant 
at  a  =  0.05? 

Goodman  and 
Kruskal’s  r 

System  Size 

0.5000 

No 

0.005 

Sample  Size 

0.2166 

No 

0.043 

Reliability 

Levels 

0.5000 

No 

0.005 

Component  Prior 
Info  Strength 

0.5000 

No 

0.005 

Component  Prior 
Info  Accuracy 

0.0004 

Yes 

0.391 

Table  7.  Pair-Wise  Factor  Interactions;  Point  Estimators. 


Factor  1 

Factor  2 

P- Value 

Significant 
at  a  =  0.05? 

Goodman  and 
Kruskal’s  r 

System  Size 

Sample  Size 

0.6368 

No 

0.053 

System  Size 

Reliability 

Levels 

0.9268 

No 

0.014 

System  Size 

Component  Prior 
Info  Strength 

0.9268 

No 

0.014 

System  Size 

Component  Prior 
Info  Accuracy 

0.0050 

Yes 

0.401 

Sample  Size 

Reliability 

Levels 

0.4014 

No 

0.092 

Sample  Size 

Component  Prior 
Info  Strength 

0.4014 

No 

0.092 

Sample  Size 

Component  Prior 
Info  Accuracy 

0.0016 

Yes 

0.478 

Reliability 

Levels 

Component  Prior 
Info  Strength 

0.9268 

No 

0.014 

Reliability 

Levels 

Component  Prior 
Info  Accuracy 

0.0050 

Yes 

0.401 

Component  Prior 
Info  Strength 

Component  Prior 
Info  Accuracy 

0.0050 

_ 

Yes 

0.401 
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Table  6  shows  that  prior  information  accuracy  is  the  only  factor  with  a  sig¬ 
nificant  main  effect.  However,  all  pair-wise  interactions  involving  prior  information 
accuracy  were  significant  as  well  (see  Table  7).  Since  the  association  between  prior 
accuracy  and  estimator  choice  (as  measured  by  Goodman  and  Kruskal’s  r)  is  an 
order  of  magnitude  larger  than  any  of  the  other  association  levels,  the  significance 
of  the  pair-wise  interactions  is  most  likely  an  artifact  of  the  dominance  of  the  prior 
accuracy  main  effect.  According  to  these  tests,  then,  prior  information  accuracy  is 
the  only  structural  characteristic  of  the  hypothetical  systems  that  affects  the  choice 
of  the  best  estimator. 

According  to  Table  5,  systems  with  inaccurate  prior  information  (the  odd- 
numbered  systems)  overwhelmingly  favored  the  MLE,  while  those  with  accurate 
component  priors  were  evenly  split  between  the  estimators.  These  results  support 
the  intuition  that,  if  there  is  a  large  gap  between  the  prior  distribution  mean  and 
the  true  component  reliability,  the  MLE  should  outperform  the  Bayes  estimator. 
However,  the  results  are  not  consistent  with  Martz  and  Waller’s  rule  of  thumb  that 
the  Bayes  estimator  outperforms  the  MLE  if  the  prior  information  is  accurate.  The 
inconsistency  may  result  from  the  fact  that,  in  this  study,  a  component’s  prior  in¬ 
formation  was  deemed  “accurate”  if  the  prior  mean  was  within  ±20  percent  of  the 
truth.  For  some  types  of  systems,  this  latitude  in  prior  information  accuracy  may 
be  too  broad. 

Results  and  Analysis:  Interval  Estimators 

Tables  8  through  11  show  the  interval  estimator  performance  measures  by  sys¬ 
tem.  The  estimator  judged  the  best  performer  for  a  given  system  is  designated  with 
an  asterisk  (*).  For  most  systems,  either  the  beta  approximation  or  the  Lindstrom- 
Madden  method  provided  the  most  desirable  interval  estimates,  although  for  three 
systems  (14,  26  and  32)  there  was  no  clearly  superior  estimator.  Appendix  C  contains 
a  discussion  of  the  Monte  Carlo  sampling  statistics  for  the  interval  estimators. 
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Table  8.  Interval  Estimator  Performance  Measures,  Systems  1-8. 


90%  Lower  Bound 


95%  Lower  Bound 


System 

R 

Method 

Avg  Limit 

Coverage 

Avg  Limit 

Coverage 

MMLI 

** 

** 

** 

L-M* 

0.5579 

1.0000 

0.4689 

1.0000 

1 

0.9960 

Beta  Approx 

0.1435 

1.0000 

0.1228 

1.0000 

BMC 

0.0952 

1.0000 

0.0784 

1.0000 

MMLI 

*♦ 

** 

** 

** 

L-M 

0.5176 

1.0000 

0.4316 

1.0000 

2 

0.9510 

0.5769 

1.0000 

0.5333 

1.0000 

BMC 

0.5114 

1.0000 

0.4609 

1.0000 

MMLI 

** 

L-M* 

0.3603 

1.0000 

0.2895 

1.0000 

Beta  Approx 

0.1932 

1.0000 

0.1752 

1.0000 

BMC 

0.1367 

1.0000 

0.1208 

1.0000 

0.5905 


MMLI 

** 

** 

L-M 

0.2256 

1.0000 

0.1726 

1.0000 

0.3686 

1.0000 

0.3424 

1.0000 

BMC 

0.2979 

1.0000 

0.2709 

1.0000 

MMLI 

** 

** 

** 

** 

L-M* 

0.4959 

1.0000 

0.4116 

1.0000 

Beta  Approx 

0.1396 

1.0000 

0.1193 

1.0000 

BMC 

0.0927 

1.0000 

0.0763 

1.0000 

0.7237 


0.7931 


MMLI 

** 

** 

** 

L-M 

0.3230 

1.0000 

0.2567 

1.0000 

0.4507 

1.0000 

0.4108 

1.0000 

BMC 

0.4029 

1.0000 

0.3600 

1.0000 

MMLI 

** 

** 

L-M-^ 

0.3947 

1.0000 

0.3202 

1.0000 

Beta  Approx 

0.1857 

1.0000 

0.1681 

1.0000 

BMC 

0.1329 

1.0000 

0.1175 

1.0000 

MMLI 

** 

** 

L-M 

0.3789 

1.0000 

0.3601 

1.0000 

0.5634 

1.0000 

0.5326 

1.0000 

BMC 

0.5357 

1.0000 

0.5027 

1.0000 

Best  Estimator. 


Not  reported. 
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Table  9.  Interval  Estimator  Performance  Measures,  Systems  9-16. 


90%  Lower  Bound 


95%  Lower  Bound 


System 

R 

Method 

Avg  Limit 

Coverage 

Avg  Limit 

Coverage 

MMLI 

0.3265 

1.0000 

0.2832 

1.0000 

L^M^^ 

0.3528 

0.9664 

0.3044 

0.9664 

9 

0.5905 

Beta  Approx 

0.1306 

1.0000 

0.1142 

1.0000 

BMC 

0.0857 

1.0000 

0.0726 

1.0000 

MMLI 

0.4515 

1.0000 

0.4006 

1.0000 

L-M 

0.5312 

0.9236 

0.4762 

1.0000 

10 

0.7738 

Beta  Approx* 

0.5849 

1.0000 

0.5482 

1.0000 

BMC 

0.5217 

1.0000 

0.4795 

1.0000 

0.9510 


MMLI 

L-M* 

0.7328 

1.0000 

0.6780 

1.0000 

Beta  Approx 

0.2669 

1.0000 

0.2469 

1.0000 

BMC 

0.2020 

1.0000 

0.1828 

1.0000 

0.9960 


MMLI 

** 

*!|C 

L-M 

0.7879 

1.0000 

0.7346 

1.0000 

0.7995 

1.0000 

0.7740 

1.0000 

BMC 

0.7611 

1.0000 

0.7292 

1.0000 

0.9684 


MMLI 

** 

** 

** 

L-M* 

0.7563 

1.0000 

0.7021 

1.0000 

Beta  Approx 

0.2581 

1.0000 

0.2324 

1.0000 

BMC 

0.1940 

1.0000 

0.1697 

1.0000 

0.8757 


MMLI 

** 

** 

** 

** 

L-M 

0.6403 

1.0000 

0.5844 

1.0000 

0.6293 

1.0000 

0.5924 

1.0000 

BMC 

0.5767 

1.0000 

0.5337 

1.0000 

MMLI 

0.3709 

1.0000 

0.3244 

1.0000 

L-M* 

0.4154 

0.9156 

0.3639 

0.9832 

Beta  Approx 

0.2067 

1.0000 

0.1895 

1.0000 

BMC 

0.1547 

1.0000 

0.1390 

1.0000 

MMLI 

** 

♦  ♦ 

** 

** 

L-M 

0.6083 

1.0000 

0.5525 

1.0000 

0.6512 

1.0000 

0.6231 

1.0000 

BMC 

0.6236 

1.0000 

0.5926 

1.0000 

0.6579 


0.8448 


Best  Estimator. 


Not  reported. 
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Table  10.  Interval  Estimator  Performance  Measures,  Systems  17-24. 


90%  Lower  Bound 


95%  Lower  Bound 


System 

R 

Method 

Avg  Limit 

Coverage 

Coverage 

MMLI 

** 

j 

** 

** 

TM^^ 

0.4721 

1.0000 

0.3899 

1.0000 

17 

0.9044 

Beta  Approx 

0.0252 

1.0000 

0.0200 

1.0000 

BMC 

0.0177 

1.0000 

0.0140 

1.0000 

MMLI 

** 

** 

** 

** 

L-M* 

0.5526 

1.0000 

0.4638 

1.0000 

18 

0.9900 

Beta  Approx 

0.3864 

1.0000 

0.3491 

1.0000 

BMC 

0.3479 

1.0000 

0.3092 

1.0000 

MMLI 

0.0654 

1.0000 

0.0448 

1.0000 

L^M^^ 

0.0856 

0.9848 

0.0593 

0.9848 

Beta  Approx 

0.0248 

1.0000 

0.0212 

1.0000 

BMC 

0.0171 

1.0000 

0.0145 

1.0000 

MMLI 

** 

** 

L-M 

0.2315 

1.0000 

0.1778 

1.0000 

0.3832 

1.0000 

0.3571 

1.0000 

1  BMC 

0.3456 

1.0000 

0.3184 

1.0000 

MMLI 

** 

** 

** 

L-M* 

0.3129 

1.0000 

0.2478 

1.0000 

Beta  Approx 

0.0184 

1.0000 

0.0144 

1.0000 

BMC 

0.0132 

1.0000 

0.0104 

1.0000 

0.5293 


MMLI 

0.1114 

1.0000 

0.0796 

1.0000 

L-M 

0.1728 

0.9320 

0.1288 

1.0000 

0.2022 

1.0000 

0.1775 

1.0000 

BMC 

0.1835 

1.0000 

0.1597 

1.0000 

MMLI 

** 

** 

** 

** 

L-M* 

0.5203 

1.0000 

0.4340 

1.0000 

Beta  Approx 

0.0481 

1.0000 

0.0420 

1.0000 

BMC 

0.0354 

1.0000 

0.0305 

1.0000 

MMLI 

** 

** 

L-M 

0.1833 

0.9262 

0.1373 

1.0000 

0.2659 

1.0000 

0.2444 

1.0000 

BMC 

0.2546 

1.0000 

0.2333 

1.0000 

Best  Estimator. 


Not  reported. 
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Table  11.  Interval  Estimator  Performance  Measures,  Systems  25-32. 


System 

R 

Method 

90%  Lower  Bound 

95%  Lower  Bound 

Avg  Limit 

Coverage 

Avg  Limit 

Coverage 

25 

0.5987 

MMLI 

0.3374 

1.0000 

0.2937 

1.0000 

L^M^ 

0.3616 

0.9624 

0.3126 

0.9624 

Beta  Approx 

0.0597 

1.0000 

0.0510 

1.0000 

BMC 

0.0449 

1.0000 

0.0379 

1.0000 

26 

0.3487 

MMLI 

0.1704 

0.9772 

0.1423 

0.9936 

L-M 

0.1546 

0.9770 

0.1227 

0.9842 

Beta  Approx 

0.1688 

1.0000 

0.1503 

1.0000 

BMC 

0.1403 

1.0000 

0.1232 

1.0000 

27 

0.9900 

MMLI 

** 

** 

** 

L-M* 

0.7820 

1.0000 

0.7284 

1.0000 

Beta  Approx 

0.0845 

1.0000 

0.0759 

1.0000 

BMC 

0.0655 

1.0000 

0.0581 

1.0000 

1 

28 

0.9044 

MMLI 

** 

** 

L-M 

0.6355 

1.0000 

0.5797 

1.0000 

l!lf4liPAl»bW«V!U 

0.6433 

1.0000 

0.6156 

1.0000 

BMC 

0.6138 

1.0000 

0.5831 

1.0000 

29 

0.6347 

MMLI 

0.3555 

1.0000 

0.3101 

1.0000 

L-M* 

0.3939 

0.9378 

0.3433 

0.9884 

Beta  Approx 

0.0468 

1.0000 

0.0396 

1.0000 

BMC 

0.0366 

1.0000 

0.0308 

1.0000 

30 

1 

0.5920 

MMLI 

0.3298 

1.0000 

0.2864 

1.0000 

L-M* 

0.3544 

0.9660 

0.3058 

0.9660 

Beta  Approx 

0.3226 

1.0000 

0.2643 

MBIIIIIIM 

BMC 

0.3012 

1.0000 

0.2725 

1.0000 

31 

0.6004 

MMLI 

0.3856 

1.0000 

0.3382 

1.0000 

L-M* 

0.4343 

0.9846 

0.3819 

0.9846 

0.0656 

1.0000 

0.0585 

1.0000 

BMC 

0.0509 

1.0000 

0.0449 

1.0000 

i 

! 

32 

0.8211 

MMLI 

** 

** 

L-M 

0.5820 

1.0000 

0.5263 

1.0000 

Beta  Approx 

0.5696 

1.0000 

0.5421 

1.0000 

BMC 

0.5466 

1.0000 

0.5172 

1.0000 

*  Best  Estimator.  **  Not  reported. 
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For  the  interval  estimators,  a  contingency  table  analysis  was  conducted  using 
methodology  similar  to  that  of  the  analysis  of  point  estimator  performance  (see  the 
previous  section);  Appendix  D  contains  the  related  contingency  tables.  As  can  be 
seen  from  Tables  12  and  13,  this  set  of  tests  yielded  the  same  conclusion  as  the  point 
estimator  tests:  only  the  accuracy  of  the  prior  information  significantly  affected  the 
choice  of  the  best  estimator  for  the  hypothetical  systems. 

One  obvious  feature  of  the  results  is  the  conservatism  of  the  interval  estimates. 
To  capture  the  degree  of  conservatism  for  each  method,  each  system’s  90  percent 
Lindstrom-Madden  and  approximate  beta  interval  estimates  were  arbitrarily  reduced 
by  10,  20  and  50  percent.  Table  14  shows  the  effect  of  the  reductions  on  the  coverage 
proportions.  The  reduction  data  show  that,  for  many  systems,  both  estimators 
produced  interval  estimates  that  were  so  conservative  they  all  would  have  covered  R 
even  if  they  were  cut  in  half. 

System  10  was  chosen  at  random  to  illustrate  the  effect  of  the  interval  reduc¬ 
tion  on  best  estimator  selection.  For  each  interval  estimation  method,  the  length 
of  each  of  System  lO’s  90  percent  interval  estimates  was  reduced  until  the  coverage 
proportions  were  equal  and  as  close  to  0.9  as  possible.  A  comparison  of  the  average 
lower  bounds  of  the  reduced  intervals  (see  Table  15)  shows  that  the  beta  approxima¬ 
tion  still  yields  the  shortest  intervals;  this  outcome  is  unchanged  from  the  original 
comparison. 
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Table  12.  Main  Factor  Effects;  Interval  Estimators. 


Factor 

P-Value 

Signihcant 
at  a  =  0.05? 

Goodman  and 
KruskaPs  r 

System  Size 

0.3612 

No 

0.016 

Sample  Size 

0.6388 

No 

0.000 

Reliability 

Levels 

0.6388 

No 

0.000 

Component  Prior 
Info  Strength 

0.3612 

No 

0.016 

Component  Prior 
Info  Accuracy 

0.0000 

Yes 

0.778 

Table  13.  Pair-Wise  Factor  Interactions;  Interval  Estimators. 


Factor  1 

Factor  2 

P- Value 

Significant 
at  a  =  0.05? 

Goodman  and 
Kruskal’s  t 

System  Size 

Sample  Size 

0.9171 

No 

0.016 

System  Size 

Reliability 

Levels 

0.9171 

No 

0.016 

System  Size 

Component  Prior 
Info  Strength 

0.6768 

No 

0.048 

System  Size 

Component  Prior 
Info  Accuracy 

0.0000 

Yes 

0.810 

Sample  Size 

Reliability 

Levels 

0.9171 

No 

0.016 

Sample  Size 

Component  Prior 
Info  Strength 

0.9171 

No 

0.016 

Sample  Size 

Component  Prior 
Info  Accuracy 

0.0000 

Yes 

0.778 

Reliability 

Levels 

Component  Prior 
Info  Strength 

0.9171 

No 

0.016 

Reliability 

Levels 

Component  Prior 
Info  Accuracy 

0.0000 

Yes 

0.778 

Component  Prior 
Info  Strength 

Component  Prior 
Info  Accuracy 

0.0000 

_ 

Yes 

0.810 
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Table  14.  Coverage  Proportions  for  Reduced  Intervals. 


System 

Method 

Reduction  Percentage 

0 

10 

30 

. 

50 

1 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

2 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

3 

L-M 

1.0000 

1.0000 

1.0000 

0.6410 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

4 

L-M 

0.8806 

0.3624 

Beta  Approx 

1.0000 

1.0000 

1.0000 

0.0000 

5 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

6 

L-M 

1.0000 

1.0000 

1.0000 

0.7218 

Beta  Approx 

1.0000 

-  _  .J 

1.0000 

1.0000 

0.3528 

~  "  n 

7 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

8 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

0.5986 

9 

L-M 

0.9664 

0.8848 

0.7474 

0.0776 

Beta  Approx 

1.0000 

1.0000 

1.0000 

0.9948 

10 

L-M 

0.9236 

0.9236 

0.9236 

0.4562 

Beta  Approx 

1.0000 

1.0000 

0.9236 

0.2404 

11 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

12 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

13 

L-M 

1.0000 

1.0000 

1.0000 

1.0000 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

14 

L-M 

1.0000 

1.0000 

1.0000 

0.7474 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

15 

L-M 

0.9156 

0.9156 

0.7724 

0.2126 

Beta  Approx 

1.0000 

1.0000 

1.0000 

1.0000 

16 

L-M 

1.0000 

1.0000 

0.8116 

0.8116 

Beta  Approx 

1.0000 

1.0000 

1.0000 

0.8116 

Table  15.  Interval  Reduction  Data  for  System  10. 


Method 

Original 

Limit 

Reduction 

Multiplier 

New 

Coverage 

New 

Avg  Limit 

Beta  Approx 

0.5849 

0.3706 

0.9236 

0.7387 

L-M 

0.5312 

0.3284 

0.9236 

0.6852 

BMC 

0.5217 

0.4395 

0.9236 

0.7320 

Conclusions 

Based  on  the  assumptions  made  and  the  results  of  the  Monte  Carlo  comparison 
study,  four  conclusions  can  be  drawn; 

1.  If  the  accuracy  of  the  component  prior  information  is  in  doubt,  then  the  MLE 
and  the  Lindstrom-Madden  estimator  are  clearly  the  best  estimators  for  R. 
However,  if  confidence  in  the  prior  information  is  high,  the  Bayes  estimators 
constructed  from  the  approximate  beta  posterior  should  be  used  instead. 

2.  Because  of  its  extreme  conservatism,  the  interval  estimator  constructed  from 
the  Bayes  Monte  Carlo  empirical  posterior  distribution  should  not  be  used  for 
component  sample  sizes  under  ten. 

3.  Easterling’s  MMLI  method  should  not  be  used  when  the  sample  sizes  are  uni¬ 
formly  under  ten,  especially  when  the  component  reliabilities  are  known  to  be 
high. 

4.  Any  small-sample  interval  estimate  generated  with  any  of  the  methods  consid¬ 
ered  will  be  conservative. 

In  summary,  the  results  of  this  investigation  show  that  both  the  beta  ap¬ 
proximation  and  the  Lindstrom-Madden  method  to  be  useful  with  small  data  sets. 
However,  neither  of  these  estimators  is  strongly  preferable  to  the  other.  ACC  could 
use  either  approach,  depending  on  the  quality  of  the  prior  information  available. 
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Recommendations  for  Further  Research 


Since  this  study  is  only  a  first  step  in  understanding  small-sample  reliability 
estimation,  ample  opportunity  exists  to  broaden  the  scope  of  the  results.  Several 
areas  with  potential  for  additional  research  are  suggested: 

1.  Additional  hypothetical  systems  should  be  considered  to  incorporate  the  fol¬ 
lowing  extensions  to  the  design: 

(a)  More  than  two  levels  of  k  and  n,  should  be  examined.  For  example,  a 
possible  extension  could  be  A:  =  {3,  5,  10,  15}  and  =  {4,  10,  20}. 

(b)  The  accuracy  and  strength  of  the  component  prior  information  and  the 
component  sample  size  should  be  allowed  to  vary  across  components. 

(c)  Parallel  and  complex  systems  should  be  considered  as  well  as  series  sys¬ 
tems. 

2.  As  mentioned  previously,  the  definition  of  “accuracy”  applied  to  the  component 
prior  information  in  this  study  may  have  been  both  simplistic  and  overly  liberal. 
The  effect  on  the  study  results  of  variations  in  prior  density  accuracy  should 
be  analyzed  more  thoroughly.  Possible  avenues  for  incorporating  error  into 
the  prior  include  adjustments  to  its  mean  or  variance  (or  both).  Such  a  study 
would  yield  valuable  information  about  the  impact  of  prior  information  on 
smaJl-sample  Bayes  estimation. 

3.  It  was  shown  earlier  in  this  chapter  that  large  decreases  in  small-sample  interval 
estimates  can  sometimes  be  made  at  a  relatively  low  cost  in  significance.  An 
analytical  understanding  of  this  phenomenon  could  result  in  the  development  of 
interval  reduction  methods,  which  would  in  turn  lead  to  small-sample  interval 
estimates  of  higher  quality. 


Appendix  A.  Bayesian  Statistical  Inference 


Introduction 

This  appendix  provides  a  brief  introduction  to  Bayesian  statistical  inference. 
With  the  aid  of  this  presentation,  the  reader  with  little  or  no  exposure  to  the  Bayesian 
paradigm  can  understand  the  Bayesian  reliability  methods  discussed  in  this  thesis. 

The  Nature  of  Probability:  Two  Views 

Before  undertaking  a  meaningful  discussion  of  Bayesian  statistical  methods, 
it  is  important  to  understand  how  and  why  those  methods  differ  from  the  classical 
statistical  techniques  taught  in  most  universities.  The  two  approaches  are  dissimilar 
because  they  are  rooted  in  different  understandings  of  the  fundamental  nature  of 
probability. 

The  classical  school  of  probability  (also  called  the  frequentist  or  objectivist 
school)  holds  that  probability  is  a  measure  of  the  relative  frequency  of  an  event 
during  repeated  experimentation.  Consider  an  experiment  for  which  Q  is  the  set  of 
all  possible  outcomes.  Define  as  one  of  those  outcomes  the  event  E  such  that  E  E  il. 
A  frequentist  statistician  would  view  the  probability  of  E  (that  is,  the  likelihood 
of  event  E  occurring)  as  P[E]  =  limn-,oo(e/n),  where  e  is  the  number  of  times  E 
occurs  during  n  experimental  trials.  Thus,  P\E]  is  the  long-run  frequency  with 
which  E  occurs.  According  to  the  classical  viewpoint,  the  value  of  this  frequency  is 
an  unknown  constant  that  can  be  reasonably  approximated  if  enough  experiments 
are  conducted. 

In  contrast,  the  Bayesian  school  (which  takes  its  name  from  the  British  clergy¬ 
man  and  amateur  mathematician  Thomeis  Bayes)  adopts  a  subjectivistic  (or  person- 
alistic)  view  of  probability.  The  Bayesian  views  P[E]  not  as  an  objective  estimate 
of  the  unknown  constant  likelihood  of  E,  but  as  a  rational  expression  of  his  or  her 
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degree  of  belief  that  E  will  occur.  Since  probability  is  subjective,  different  people 
can  assign  different  values  to  P[E]  based  on  varying  perceptions  of  prior  experience. 
Thus,  from  the  Bayesian  perspective,  this  prior  experience  replaces  experimentation 
as  the  source  of  the  probability  estimate.  However,  the  Bayesian  does  not  abandon 
experimentation,  because  it  is  an  important  means  of  refining  his  or  her  assessment 
of  P[E]. 

To  illustrate  these  divergent  views  of  the  nature  of  probability,  consider  a 
simple  experiment.  A  certain  coin  is  flipped,  with  the  outcome  being  either  “heaxis” 
or  “tails”;  that  is,  0  =  {heads,  tails}.  Suppose  the  likelihood  of  each  outcome  is 
unknown.  Further,  define  the  event  E  =  {heads}.  In  order  to  make  a  probability 
statement  about  E,  the  frequentist  would  perform  n  experiments  (where  n  is  large, 
say  1000),  observe  m  occurrences  of  E,  and  declare  that,  for  this  particular  coin, 
P[E]  «  m/n.  On  the  other  hand,  the  Bayesian  might  conclude,  based  on  experience 
with  similar  coins,  that  E  is  just  as  likely  to  occur  as  not.  Therefore,  he  or  she  would 
declare  that  P[E]  =  1  /2,  then  conduct  a  series  of  experimental  trials  to  validate  and 
refine  that  estimate. 

There  is  a  great  deal  of  controversy  over  the  legitimacy  of  the  Bayesi2m  view 
of  probability.  Many  theoretical  statisticians  reject  the  Bayesian  approach  because 
they  believe  that  statistics,  as  a  science,  should  be  as  objective  as  possible.  However, 
Press  and  others  have  pointed  out  that  the  methods  of  classical  statistical  inference 
incorporate  a  great  deal  of  subjectivity,  although  they  purport  to  be  objective.  The 
subjective  nature  of  the  Bayesian  statistics  is  not  a  weakness  but  a  strength,  because 
valuable  p^ior  experience  can  be  accounted  for  in  the  inference  process  [36:1-14]. 

Bayes’  Theorem 

Bayes’  Theorem,  which  is  presented  below,  is  both  stated  and  proven  using 
the  definition  of  conditional  probability.  This  definition,  which  can  be  found  in  any 
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introductory  probability  text,  provides  a  convenient  mechanism  for  incorporating 
subjective  probability  information  into  a  statistical  analysis. 


Theorem  A.l  (Bayes)  Let  ft  be  a  sample  space.  Define  the  mutually  exclusive 
events  Fi  such  that  (1)  Fi  6  ft  Vt  and  (2)  U^i  Further,  let  E  be  an 

event  such  that  (1)  the  conditional  probabilities  P[E\Fi]  are  defined  for  all  i,  and  (2) 
P[E]  >  0.  Then 


E”,  f 


The  proof  of  this  theorem,  which  is  a  straightforward  application  of  the  definition  of 
conditional  probability,  is  given  by  Biswas  [5:28]. 

Bayes’  Theorem  had  its  origin  in  Bayes’  1763  paper  An  Essay  Towards  Solving 
a  Problem  in  the  Doctrine  of  Chances,  where  it  was  originally  limited  to  the  binomial 
sampling  problem  [3].  Laplace,  working  independently,  published  the  general  form 
of  this  theorem  in  1774  [18];  he  developed  it  more  ful  y  in  several  subsequent  papers. 
The  modern  fusion  of  Bayes’  Theorem  with  subjective  probability  is  due  mostly  to 
the  1939  work  of  Jeffreys;  De  Groot,  Savage,  de  Finetti,  Ramsey  and  Renyi  also 
made  significant  contributions  to  the  development  of  Bayesian  theory  [36:8,16]. 


Foundations  of  Bayesian  Inference 

Consider  a  failure  density  with  some  parameter  9.  Based  on  earlier  discussion, 
the  frequentist  view  of  9  is  that  its  true  value  is  a  fixed,  unknown  constant  that  can 
be  discovered  through  infinite  sampling.  From  the  Bayesian  perspective,  however, 
the  true  value  of  9  could  be  any  of  a  range  of  values,  some  of  which  are  more  likely 
to  be  the  true  value  of  9  than  others.  It  is  not  surprising,  then,  that  the  perceived 
value  of  9  is  viewed  as  a  random  variable  in  Bayesian  inference.  (In  the  discussion 
that  follows,  the  symbol  0  will  be  used  to  refer  to  the  random  variable;  9  will  refer 
to  a  particular  value  of  0.) 
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Since  0  is  viewed  as  a  random  variable,  subjective  information  about  the  likeli¬ 
hood  of  its  possible  values  may  be  encoded  in  a  probability  density  or  mass  function 
g{0).  Once  g{6)  is  defined,  experiments  can  be  conducted  to  validate  and  refine  the 
assessment  of  0.  Bayes’  Theorem  is  the  mechanism  for  combining  the  subjective 
information  with  the  experimental  results. 

The  subjective  assessment  of  0  is  refined  in  a  straightforward  manner.  First, 
a  random  sample  of  size  n  of  independent,  identically-distributed  random  variables 
Xi  is  taken  ({X,  =  0:^,1  =  1 . .  .n}).  The  probability  mass/density  function  of  X,  is 
f{xi\9)  (that  is,  the  form  of  the  density /mass  function  of  X,  depends  on  the  value  of 
0).  After  the  random  sample  is  taken,  the  information  from  the  sample  is  combined 
with  the  a  priori  information  assumed  about  0  (i.  e.,  g{9))  to  get  an  d  posteriori 
refined  probability  statement  about  0.  If  0  is  assumed  to  be  a  discrete  random 
variable,  this  procedure  is  captured  in  the  following  extension  to  Bayes’  Theorem. 


Theorem  A. 2  Let  Q  be  a  discrete  random  variable  with  probability  mass  function 
g{6).  Further,  let  X  =  {Xi,...,Xn}  be  a  set  of  independent,  identically-distributed 
random  variables,  each  with  probability  mass/density  function  f(x\0),  and  let  x  = 
{xi, . . . ,  x„}  be  a  set  of  random  samples  from  those  random  variables.  Then 


h{0\x) 


The  proof  of  this  theorem  is  omitted  because  it  is  trivial  [36:24].  A  similar  extension 
to  Bayes’  Theorem  can  be  derived  for  continuous  0.  (Note  that  X;  can  be  either 
discrete  or  continuous,  despite  of  the  form  of  g{0).) 

Since  g{0)  incorporates  d  priori  information  about  0,  it  is  called  the  prior 
density  of  0;  its  parameters  are  called  hyperparameters  to  differentiate  them  from 
the  parameters  of,  say,  f{x\6).  Similarly,  h{0\x)  is  called  the  posterior  density  of 
0.  These  respective  terms  are  usually  shortened  to  “prior”  and  “posterior”  for 
simplicity. 
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The  Prior  Distribution 


It  should  be  obvious  that  a  Bayesian  estimate  is  only  aa  good  as  the  accuracy 
of  the  prior  assessment  of  0.  This  means,  then,  that  the  selection  of  the  prior  density 
g{6)  must  be  made  with  great  care. 

The  prior  density  of  0  can  take  any  form  that  seems  reasonable.  Usually,  a 
proper  prior  (a  probability  density  or  mass  function)  is  used.  However,  the  prior  need 
not  integrate  to  unity  (in  which  case  it  is  called  an  improper  prior).  The  prior  can 
even  be  an  empirical  function.  In  most  practical  applications,  however,  priors  are 
selected  from  two  classes:  natural  conjugate  priors  (sometimes  called  convenience 
priors),  and  noninformative  priors  (also  called  indifference  or  vague  priors). 

When  a  natural  conjugate  prior  is  used,  the  resulting  posterior  is  of  the  same 
family.  The  form  of  this  prior  depends  on  the  form  of  f{x\6).  For  example,  a  natural 
conjugate  prior  of  a  binomial  probability  is  the  beta  density.  This  prior  returns 
a  posterior  that  is  also  a  beta  density.  Press  discusses  the  advantages  of  natural 
conjugate  priors  at  length,  and  catalogs  the  common  conjugate  prior  families  [36]. 

A  noninformative  prior  is  used  when  there  is  little  or  no  reason  to  prefer  any  one 
value  of  0  over  another.  There  are  various  types  of  noninformative  priors  discussed 
in  the  literature;  the  most  common  form  used  is  the  uniform  distribution.  The  use 
of  noninformative  priors  has  been  the  subject  of  much  debate;  many  statisticians  see 
no  justification  for  the  Bayesian  approach  if  no  meaningful  prior  information  exists. 
See  Press  [36:46-51]  for  a  detailed  discussion  of  the  noninformative  prior. 

Once  the  family  of  the  prior  is  chosen,  its  hyperparameters  can  be  assigned 
values  based  on  any  reasonable  criterion.  In  some  cases  it  is  perfectly  acceptable  to 
choose  hyperparameter  values  based  on  “gut  feel.”  For  example,  a  group  of  design 
engineers  might  feel  strongly  (perhaps  because  of  prior  design  or  manufacturing  ex¬ 
perience)  that  a  certain  switch  has  a  failure  probability  of  0.01.  In  other  situations, 
it  may  be  advisable  to  estimate  the  hyperparameters  from  a  data  set  derived  from 
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a  different  (but  conceptually  related)  experimental  situation.  In  any  case,  the  hy¬ 
perparameters  should  be  estimated  with  care.  If  the  choice  is  made  blindly,  it  can 
render  the  entire  analysis  meaningless.  For  example,  Kapur  and  Lamberson  show 
that  the  parameters  of  a  beta  prior  can  be  selected  so  that  the  experimental  results 
will  not  substantially  affect  the  form  of  the  posterior,  even  if  those  results  completely 
contradict  the  the  prior  information  [15:380-382]. 

Determining  the  form  of  g{6)  is  a  difficult  process,  and  the  results  are  often 
controversial.  The  key  point  to  remember  when  selecting  the  prior  is  that,  although 
the  choice  is  ultimately  based  on  subjective  factors,  that  choice  should  be  both 
reasonable  and  justifiable. 

Bayesian  Point  Estimation 

A  point  estimate  is  essentially  a  decision  about  the  value  of  the  parameter  being 
estimated.  It  is  possible,  if  not  rigorous,  to  avoid  decision  theory  when  introducing 
the  basics  of  classical  point  estimation.  However,  Bayesian  point  estimation  must  be 
introduced  in  a  decision-theoretic  framework. 

Before  Bayesian  point  estimation  is  discussed,  several  concepts  from  decision 
theory  must  be  defined.  These  definitions,  which  deal  with  a  single  parameter  only, 
are  adapted  from  the  general  forms  given  by  Mood  et  al.  [33:297-299,  344]. 

Definition  A.l  Let  0  =  0(Xi,...,X„)  be  an  estimator  of  a  parameter  6,  and 
let  9  be  a  particular  value  of  0.  The  loss  function  measures  the  error  sustained 
by  estimating  6  by  0.  This  function  is  denoted  by  X{9',0),  and  has  the  following 
properties:  (1)  \{9]6)  is  a  real  function,  (2)  X{9',9)  >  0  and  (3)  X{9',d)  =  0 
when  0  =  9. 
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There  are  several  possible  loss  functions.  The  one  most  often  seen  in  practical 
applications  is  the  squared-error  loss  function: 

\{e-e)  =  (e-of  (46) 

Definition  A. 2  For  some  loss  function  0),  the  risk  function  of  the  estimator  0 
is 

Pe(^)  =  £lA(0;^)] 

The  risk  associated  with  a  squared-error  loss  function  is  E[(0  ~  which  is  easily 
recognized  as  the  mean  square  error  of  0. 

Definition  A. 3  The  Bayes  risk  of  the  estimator  0  with  respect  to  the  loss  function 
\{0]  0)  and  the  prior  density  g{0)  is 

^A..(e)  =  /  P^i0)gi0)d0 

«/© 

With  the  above  terms  defined,  Ba"es  estimators  can  now  be  introduced. 

Definition  A. 4  The  Bayes  estimator  0*  of  0  (with  respect  to  the  loss  function 
A(0;0)  and  the  prior  density  g{0))  is  the  estimator  with  the  smallest  Bayes  risk;  that 
is,  0*  is  the  estimator  for  which 

Ra,j(0‘)  =  inf(1JA,,(e)} 
e 

It  is  straightforward  to  find  the  Bayes  estimator  for  a  given  loss  function  and 
prior  density.  For  a  squared-error  loss  function,  Mood  et  al.  show  that  the  Bayes 
estimator  is  the  mean  of  the  posterior  density  (i.e.,  0*  =  E[0|x])  [33:345]. 

One  other  property  of  Bayes  estimators  is  useful  in  certain  situations.  Biswas 
shows  that  Bayes  estimates  and  maximum  likelihood  estimates  converge  as  sample 
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size  approaches  infinity.  Further,  this  property  holds  regardless  of  the  form  of  the 
prior  [5:338]. 

Bayesian  Interval  Estimation 

The  interpretation  of  a  Bayesian  interval  estimate  is  different  from  that  of 
its  classical  counterpart,  the  confidence  interval.  RecaJl  that  a  frequentist  assumes 
the  parameter  9  estimated  is  an  unknown  constant.  In  a  classical  interval  estimation 
method,  the  estimators  of  the  end  points  of  the  interval  are  random  variables.  There¬ 
fore,  a  classical  confidence  coefficient  is  not  an  explicit  probability  statement  about 
9.  Instead,  it  is  the  probability  that  the  interval  estimator  will  generate  an  interval 
that  will  contain  the  true  value  of  9  (and  is  therefore  a  declaration  of  confidence  in 
the  estimate)  [32:353]. 

In  the  Bayesian  view,  however,  0  is  a  random  variable,  so  a  Bayesian  proba¬ 
bility  (or  credibility)  interval  is  a  direct  probability  statement  about  the  value  of  ©. 
Thus,  the  two-sided  100(1  —  7)  percent  probability  interval  \9i,9u\  has  an  intuitive 
interpretation:  0  will  take  on  a  value  within  the  interval  with  probability  1  —  7  (that 
is,  P{0  €  [9li9u\}  =  1—7).  The  bounds  of  the  interval  [9i,9u]  can  be  evaluated 
directly  by  solving 

r  =  I  (47) 

»/— 00  z 

and 

hWi)d9  =  I  (48) 

for  9i,  and  9u-  This  leads  to  a  symmetric  probability  interval  for  0.  There  are  other 
types  of  two-sided  Bayesian  probability  intervals,  but  the  symmetric  interval  is  most 
often  used  because  it  is  simple  to  calculate  [30:208]. 

A  one-sided  Bayesian  probability  interval  has  an  interpretation  similar  to  that 
of  the  two-sided  interval.  The  lower  probability  bound,  which  is  considered  exten- 
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sively  in  this  thesis,  is  found  by  solving  the  following  equation  for  Ol'. 

h{e\x)de  =  7  (49) 

J  —  QO 

This  is  an  obvious  adaptation  of  Equation  (47).  An  expression  for  the  upper  bound 
can  be  formulated  by  similarly  manipulating  Equation  (48). 

Empirical  Bayesian  Inference 

A  school  of  thought  that  tries  to  mitigate  the  subjectivity  of  Bayesian  inference 
is  the  empirical  Bayes  approach.  Empirical  Bayes  procedures  compensate  for  the  un¬ 
certainty  over  the  form  of  the  prior  density  by  estimating  its  hyperparameters  from 
the  current  data  set.  Because  they  violate  Bayes  Theorem  by  making  the  hyperpa¬ 
rameters  dependent  on  the  test  data,  empirical  Bayes  procedures  are  controversial 
[36:43].  Furthermore,  empirical  Bayes  estimates  require  a  data  set  of  adequate  size; 
when  the  data  set  is  small,  it  makes  no  sense  to  use  these  methods  (philosophical 
issues  notwithstanding).  Since  empirical  procedures  are  not  used  in  this  thesis,  they 
are  not  developed  here.  The  interested  reader  is  referred  to  [28]  and  Chapter  13  of 
(30], 

Further  Reading 

Press,  who  has  been  cited  frequently  in  this  appendix,  provides  a  general  intro¬ 
duction  to  Bayesian  statistics  that  is  well- written  and  easy  to  understand  [36].  Also, 
Martz  and  Waller’s  reliability  text  [30]  and  Berger’s  book  on  decision  theory  [4]  con¬ 
tain  brief,  but  comprehensible,  introductions  to  Bayesian  inference  as  a  foundation 
for  their  primary  subject  matter  (with  Berger’s  treatment  being  considerably  more 
rigorous).  For  an  advanced  treatment  of  Bayesian  probability  theory,  the  reader 
should  refer  to  De  Groot  [10]. 
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Most  standard  reliability  textbooks  contain  at  least  a  chapter  introducing 
Bayesian  inference  and  describing  its  applications  to  reliability.  Some  of  the  texts 
containing  particularly  enlightening  expositions  of  Bayesian  theory  include  Crowder 
et  al.  [9],  Kapur  and  Lamberson  [15],  and  Mann  et  al.  [27]. 

For  a  reliability  text  written  exclusively  from  the  Bayesian  viewpoint,  the 
reader  is  referred  to  Martz  and  Waller  [30].  This  book  is  the  standard  reference 
work  on  Bayesian  reliability  methods. 
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Appendix  B.  The  Mellin  Transform 

The  Mellin  transform  of  f{x)  with  respect  to  the  complex  parameter  u  is 

too 

M{/(x);u}  =  /  (50) 

Jo 

If  /(x)  is  the  density  function  of  the  continuous  random  variable  X  and  u  is  an 
integer  greater  than  2,  it  immediately  follows  that 

M{/(x);«}  =  £;[X“-^]  (51) 

which  is  the  {u  —  l)th  central  moment  of  X. 

The  Mellin  transform  is  useful  for  finding  density  functions  for  products  and 
quotients  of  random  variables.  Consider  the  independent,  nonnegative  random  vari¬ 
ables  Ai, . . . ,  Afc  with  density  functions  /i(xi), . . . ,  /jt(xjt).  If  the  random  variable 
Y  =  nfsi  Xi  has  a  density  function  h{y),  it  can  be  shown  that 

M{h{yy,u}  =  J|M{/,(x,);u}  (52) 

t=i 

The  function  /(x)  associated  with  the  Mellin  transform  M{/(x);«}  can  be 
derived  using  the  Mellin  inversion  integral: 

1  r  /•p+«» 

f{^)  =  ^  [jim  X  “M{/(x);u}duJ  (53) 

Since  tables  of  Mellin  transform  pairs  are  available  for  common  functional  forms,  it 
is  often  possible  to  avoid  evaluating  the  inversion  integral. 

Springer,  from  whose  book  this  discussion  is  digested,  provides  a  comprehensive 
treatment  of  the  Mellin  transform  [40:91-97].  Oberhettinger  provides  an  extensive 
catalog  of  Mellin  transform  pairs  [35]. 


Appendix  C.  Assessment  of  the  Monte  Carlo  Process 

Overview 

When  a  Monte  Carlo  study  is  conducted,  it  is  important  to  understand  the 
variability  inherent  in  the  simulation  process.  If  a  simulation’s  variance  is  too  high, 
its  results  are  suspect;  this  usually  happens  when  the  sample  size  of  the  simulation  is 
too  small.  Two  measures  of  merit  are  often  used  to  judge  the  quality  of  a  Monte  Carlo 
simulation:  the  Monte  Carlo  standard  error  (MCSE)  and  the  Monte  Carlo  confidence 
interval  (MCCI).  This  appendix  presents  the  theory  behind  these  statistics,  then  uses 
them  to  assess  the  simulation  conducted  during  this  investigation. 

Theory 

Since  the  variability  of  the  simulation  process  is  under  consideration,  the  MCSE 
is  the  standard  error  constructed  from  the  expected  value  of  the  parameter  being 
simulated.  Consider  a  simulation  parameter  9  with  mean  fig  and  variance  aj.  Fur¬ 
ther,  define  N  as  the  Monte  Carlo  sample  size.  From  Theorem  7.1  in  Mendenhall 
et  al.  [32:304],  it  follows  that  the  Monte  Carlo  variance  is 

2 

<^MC  =  Var[pg]  =  ^  (54) 

so  the  MCSE  is  obviously 

(tmc  =  (55) 

The  unbiased  estimator  for  (Tmc  is 
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where  {0,,  i  =  1, . . . ,  iV}  is  the  set  of  simulated  estimates  and 


H  = 


N 


(57) 


is  the  estimator  of  fig. 

If  the  MCSE  is  known,  the  general  form  of  the  MCCI  can  easily  be  derived 
using  the  pivot  method.  Since  fig  is  normally  distributed  with  standard  error  (Tmci 
it  immediately  follows  that  the  100(1  —  7)  percent  MCCI  for  fig  is 

[fig  ±  •  (Tmc]  (58) 

Results  and  Analysis 

Tables  16  lists  the  estimated  MCSE  and  a  95  percent  MCCI  by  system  for  each 
type  of  point  estimate.  Table  17  contains  the  same  information  for  lower  confidence 
and  probability  limits  generated  by  the  Lindstrom-Madden  and  approximate  beta 
interval  estimators  (due  to  technical  difficulties,  no  statistics  are  presented  for  the 
other  two  interval  estimators).  The  sampling  statistics  show  that,  for  N  —  5000,  the 
variability  in  the  simulation  process  is  reasonably  low. 

The  Monte  Carlo  statistics  for  the  Bayesian  estimators  are  much  lower  than 
those  derived  for  the  classical  estimators.  However,  this  is  not  surprising  since  the 
prior  information  dominates  the  small  sets  of  simulated  test  data.  A  comparison  of 
the  Bayesian  results  in  Table  16  with  the  MSE  values  for  the  Bayes  estimator  (see 
Table  5)  suggests  that  there  is  a  significant  bias  component  in  the  Bayes  estima¬ 
tor’s  MSE;  the  size  of  the  bias  seems  to  be  a  function  of  the  accuracy  of  the  prior 
information. 
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Table  16.  Monte  Carlo  Statistics:  Point  Estimates. 


System 

MLE 

Bayes  Estimator 

Smc 

95%  MCCI 

Smc 

95%  MCCI 

1 

0.0005 

[0.9945,  0.9964] 

0.0000 

[0.2369,  0.2370] 

2 

0.0015 

[0.9498,  0.9557] 

0.0003 

[0.7143,  0.7155] 

3 

0.0029 

[0.7663,  0.7776] 

0.0002 

[0.2663,  0.2669] 

4 

0.0032 

■fimintnBirai 

0.0003 

[0.4646,  0.4659] 

5 

0.0018 

0.0001 

[0.2314,  0.2320] 

6 

0.0030 

[0.7185,  0.7305] 

0.0006 

[0.5890,  0.5914] 

7 

0.0027 

[0.8078,  0.8185] 

0.0001 

[0.2573,  0.2579] 

8 

0.0028 

[0.7890,  0.8000] 

0.0003 

[0.6650,  0.6662] 

9 

0.0020 

[0.5860,  0.5939] 

0.0004 

[0.2026,  0.2041] 

10 

0.0018 

[0.7705,  0.7775] 

0.0008 

[0.7011,  0.7040] 

11 

0.0010 

[0.9498,  0.9536] 

0.0001 

[0.3445,  0.3449] 

12 

0.0003 

[0.9945,  0.9957] 

0.0001 

[0.8706,  0.8709] 

13 

0.0008 

[0.9689,  0.9719] 

0.0002 

[0.3603,  0.3609] 

14 

■miiiEk 

[0.8714,  0.8770] 

0.0006 

[0.7428,  0.7450] 

15 

0.0020 

(0.6536,  0.6614] 

0.0003 

[0.2749,  0.2759] 

16 

0.0016 

[0.8424,  0.8486] 

0.0004 

[0.7397,  0.7413] 

17 

0.0021 

[0.8994,  0.9074] 

0.0000 

[0.0544,  0.0545] 

18 

0.0007 

[0.9885,  0.9912] 

0.0001 

[0.5228,  0.5232] 

19 

0.0027 

[0.3436,  0.3544] 

0.0000 

[0.0423,  0.0425] 

20 

0.0032 

[0.5934,  0.6059] 

0.0003 

[0.4784,  0.4797] 

21 

0.0031 

[0.7055,  0.7176] 

0.0001 

[0.0420,  0.0422] 

22 

0.0032 

[0.5026,  0.5150] 

0.0006 

[0.3056,  0.3073] 

23 

0.0014 

[0.9529,  0.9586] 

0.0000 

[0.0767,  0.0768] 

24 

0.0031 

[0.5208,  0.5330] 

[0.3492,  0.3503] 

25 

0.0020 

[0.5962,  0.6039] 

0.0002 

[0.1012,  0.1019] 

26 

0.0017 

[0.3437,  0.3502] 

0.0005 

[0.2457,  0.2477] 

27 

0.0004 

[0.9897,  0.9914] 

0.0000 

[0.1224,  0.1224] 

28 

0.0015 

[0.8670,  0.8727] 

0.0003 

[0.7311,  0.7324] 

29 

0.0020 

[0.6308,  0.6387] 

0.0001 

[0.0821,  0.0827] 

30 

0.0020 

[0.5881,  0.5959] 

0.0007 

[0.4293,  0.4318] 

31 

[0.6740,  0.6815] 

[0.0973,  0.0976] 

32 

0.0017 

[0.8184,  0.8249] 

0.0004 

[0.6607,  0.6621] 
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Table  17.  Monte  Carlo  Statistics:  Interval  Lower  Bounds  (Selected  Estimators) 


System 

Lindstrom-Madden  Limit 

Approximate  Beta  Limit 

Smc 

95%  MCCI 

Smc 

95%  MCCI 

1 

0.0005 

[0.5570,  0.5589] 

0.0000 

[0.1434,  0.1435] 

2 

0.0014 

[0.5148,  0.5203] 

0.0003 

[0.5763,  0.5776] 

3 

0.0024 

[0.3556,  0.3649] 

0.0001 

[0.1930,  0.1935] 

4 

0.0023 

[0.2211,  0.2300] 

0.0003 

[0.3680,  0.3692] 

5 

0.0016 

[0.4926,  0.4991] 

0.0001 

[0.1394,  0.1398] 

6 

0.0024 

[0.3183,  0.3277] 

0.0006 

[0.4495,  0.4519] 

7 

0.0023 

[0.3901,  0.3992] 

0.0001 

[0.1855,  0.1859] 

8 

0.0023 

[0.3744,  0.3835] 

0.0003 

[0.5628,  0.5641] 

9 

0.0018 

[0.3492,  0.3565] 

0.0003 

[0.1301,  0.1312] 

10 

0.0019 

[0.5275,  0.5349] 

0.0008 

[0.5833,  0.5865] 

11 

0.0012 

[0.7305,  0.7351] 

0.0001 

[0.2667,  0.2671] 

12 

0.0004 

[0.7872,  0.7887] 

0.0001 

[0.7993,  0.7996] 

13 

0.0010 

[0.7544,  0.7582] 

0.0001 

[0.2578,  0.2583] 

14 

0.0016 

[0.6371,  0.6436] 

0.0006 

[0.6281,  0.6305] 

15 

0.0019 

0.0002 

16 

0.0018 

[0.6049,  0.6118] 

0.0004 

[0.6504,  0.6521] 

17 

0.0019 

[0.4684,  0.4758] 

0.0000 

[0.0251,  0.0252] 

18 

0.0007 

[0.5512,  0.5539] 

0.0001 

[0.3862,  0.3866] 

19 

0.0015 

[0.0827,  0.0884] 

0.0000 

[0.0247,  0.0248] 

20 

0.0023 

[0.2270,  0.2360] 

0.0003 

[0.3826,  0.3838] 

21 

0.0024 

[0.3082,  0.3177] 

0.0000 

[0.0184,  0.0185] 

22 

0.0021 

[0.1687,  0.1769] 

0.0004 

[0.2015,  0.2029] 

23 

0.0014 

[0.5176,  0.5229] 

0.0000 

[0.0481,  0.0482] 

24 

0.0021 

[0.1792,  0.1874] 

0.0002 

[0.2654,  0.2664] 

25 

0.0018 

[0.3581,  0.3652] 

0.0001 

[0.0595,  0.0600] 

26 

0.0012 

[0.1522,  0.1570] 

0.0004 

[0.1680,  0.1696] 

27 

0.0006 

[0.7809,  0.7831] 

0.0000 

[0.0846,  0.0846] 

28 

0.0017 

[0.6323,  0.6388] 

0.0004 

[0.6426,  0.6439] 

29 

0.0019 

[0.3902,  0.3976] 

0.0001 

[0.0466,  0.0470] 

30 

0.0018 

[0.3509,  0.3580] 

0.0006 

[0.3215,  0.3237] 

31 

0.0019 

[0.4306,  0.4380] 

0.0001 

[0.0655,  0.0658] 

32 

0.0018 

[0.5785,  0.5856] 

0.0004 

[0.5689,  0.5704] 
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Appendix  D.  Contingency  Tables 


Point  Estimator  Tables 

This  section  contains  the  contingency  tables  that  were  used  to  test  for  the 
significance  of  factor  effects  on  the  choice  of  the  best  point  estimator.  The  single¬ 
factor  tables  (numbers  18  to  22)  are  listed  first,  followed  by  the  pair-wise  tables 
(numbers  23  to  32).  For  two-by-two  tables,  Fisher’s  Exact  Test  was  used  to  derive 
the  p-values;  for  larger  tables,  the  p-value  was  found  using  Pearson’s  chi-squared 
statistic.  The  values  for  Goodman  and  Kruskal’s  t,  a  measure  of  nominal  association, 
are  also  included. 


Table  18.  Main  Effect  Test:  System  Size. 


MLE 

Bayes 

Row  Totals 

fc  =  5 

11 

5 

16 

o 

11 

12 

4 

16 

Column  Totals 

23 

9 

32 

P- Value  =  0.5000  t  =  0.005 


Table  19.  Main  Effect  Test:  Sample  Size. 


MLE 

Bayes 

Row  Totals 

n,  =  4 

10 

6 

16 

Tli  =  10 

13 

3 

16 

Column  Totals 

23 

9 

32 

P- Value  =  0.2166  r  =  0.043 
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Table  20.  Main  Effect  Test;  /2,  Levels. 


MLE 

Bayes 

Row  Totals 

equal 

12 

IKH 

16 

unequal 

11 

5 

16 

Column  Totals 

23 

9 

32 

P- Value  =  0.5000  r  =  0.005 


Table  21.  Main  Effect  Test:  Prior  Information  Strength. 


MLE 

Bayes 

Row  Totals 

uniformly  weak 

12 

4 

16 

uniformly  strong 

11 

5 

16 

Column  Totals 

23 

9 

32 

P- Value  =  0.5000  t  =  0.005 


Table  22.  Main  Effect  Test:  Prior  Information  Accuracy. 


MLE 

Bayes 

Row  Totals 

uniformly  inaccurate 

16 

0 

16 

uniformly  accurate 

7 

9 

16 

Column  Totals 

23 

9 

32 

P- Value  =  0.0004  r  =  0.391 
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Table  23.  Interaction  Test:  System  Size  and  Sample  Size. 


MLE 

Bayes 

Row  Totals 

5 

3 

8 

k  =  5  and  n,  =  10 

6 

2 

8 

5 

3 

8 

7 

1 

8 

Column  Totals 

23 

9 

32 

.Y2  1.700  P- Value  =  0.6368  t  =  0.053 


Table  24.  Interaction  Test:  System  Size  and  Hi  Levels. 


MLE 

Bayes 

Row  Totals 

fc  =  5  and  R,  equal 

6 

2 

8 

fc  =  5  and  R,  unequal 

5 

3 

8 

A:  =  10  and  R,  equal 

6 

2 

8 

A;  =  10  and  Ri  unequal 

6 

2 

8 

Column  Totals 

23 

' 

9 

32 

=  0.464  P- Value  =  0.9268  r  =  0.014 


Table  25.  Interaction  Test:  System  Size  and  Prior  Information  Strength. 


MLE 

Bayes 

Row  Totals 

k  =  5  and  weak  priors 

6 

2 

8 

A;  =  5  and  strong  priors 

5 

3 

8 

A:  =  10  and  weak  priors 

6 

2 

8 

A:  =  10  and  strong  priors 

6 

2 

8 

Column  Totals 

23 

9 

32 

I 


=  0.464  P- Value  =  0.9268  r  =  0.014 


Table  26.  Interaction  Test:  System  Size  and  Prior  Information  Accuracy. 


MLE 

Bayes 

Row  Totals 

k  =  5  and  inaccurate  priors 

8 

0 

8 

k  =  5  and  accurate  priors 

3 

C 

u 

8 

=  10  and  inaccurate  priors 

8 

0 

8 

k  —  10  and  accurate  priors 

4 

4 

8 

Column  Totals 

23 

9 

32 

=  12.831  P-Value  =  0.0050  t  =  0.401 


Table  27.  Interaction  Test:  Sample  Size  and  7?,  Levels. 


MLE 

Bayes 

Row  Totals 

n,  =  4  and  R,  equal 

6 

2 

8 

=  4  and  Ri  unequal 

4 

4 

8 

n,  =  10  and  R,  equal 

6 

2 

8 

Tii  =  10  and  Ri  unequal 

V 

1 

8 

Column  Totals 

23 

9 

32 

=  2.937  P-Value  =  0.4014  r  =  0.092 


Table  28.  Interaction  Test:  Sample  Size  and  Prior  Information  Strength. 


MLE 

Bayes 

Row  Totals 

n,  =  4  and  weak  priors 

6 

2 

8 

n,  =  4  and  strong  priors 

4 

4 

8 

n,  =  10  and  weak  priors 

6 

2 

8 

n,  =  10  and  strong  priors 

7 

1 

8 

Column  Totals 

23 

9 

32 

X^  =  2.937  P-Value  =  0.4014  r  =  0.092 
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Table  29,  Interaction  Test;  Sample  Size  and  Prior  Information  Accuracy. 


MLE 

Bayes 

Row  Totals 

Ui  =  4  and  inaccurate  priors 

8 

0 

8 

Tii  =  4  and  accurate  priors 

2 

6 

8 

nj  =  10  and  inaccurate  priors 

8 

0 

8 

ni  =  10  and  accurate  priors 

5 

3 

8 

Column  Totals  23 

9 

32 

=  15.304  P-Value  =  0.0016  r  =  0.478 


Table  30.  Interaction  Test:  72,  Levels  and  Prior  Information  Strength. 


MLE 

Bayes 

Row  Totals 

Ri  equal  and  weak  priors 

6 

2 

8 

Ri  equal  and  strong  priors 

6 

2 

8 

72,  unequal  and  weak  priors 

6 

2 

8 

72,  unequal  and  strong  priors 

5 

3 

8 

Column  Totals 

23 

9 

32 

=  0.464  P- Value  =  0.9268  r  =  0.014 


Table  31.  Interaction  Test:  72,  Levels  and  Prior  Information  Accuracy. 


MLE 

Bayes 

72,  equal  and  inaccurate  priors 

8 

0 

8 

72,  equal  and  accurate  priors 

4 

4 

8 

Ri  unequal  and  inaccurate  priors 

8 

0 

8 

Ri  unequal  and  accurate  priors 

3 

5 

8 

Column  Totals 

23 

9 

32 

X^  =  12.831  P- Value  =  0.0050  t  =  0.401 
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Table  32.  Interaction  Test:  Prior  Information  Strength  and  Accuracy. 


MLE 

Bayes 

Row  Totals 

weak,  inaccurate  priors 

8 

0 

8 

weak,  accurate  priors 

4 

4 

8 

strong,  inaccurate  priors 

8 

0 

8 

strong,  accurate  priors 

3 

5 

8 

Column  Totals 

23 

9 

32 

=  12.831  P- Value  =  0.0050  r  =  0.401 


Interval  Estimator  Tables 

This  section  contains  the  contingency  tables  that  were  used  to  test  for  the 
significance  of  factor  effects  on  the  choice  of  the  best  interval  estimator.  Frequency 
counts  for  the  tied  caises  (Systems  14,  26  and  32)  were  grouped  with  the  frequency 
counts  of  those  systems  preferring  the  approximate  beta  estimator;  this  was  necessary 
to  keep  the  expected  cell  frequencies  sufficiently  large. 

The  single-factor  tables  (numbers  33  to  37)  are  listed  first,  followed  by  the 
pair-wise  tables  (numbers  38  to  47).  For  two-by-two  tables,  Fisher’s  Exact  Test  was 
used  to  derive  the  p-values;  for  larger  tables,  the  p-value  was  found  using  Pearson’s 
chi-squared  statistic.  The  values  for  Goodman  and  Kruskal’s  t,  a  measure  of  nominal 
association,  are  also  included. 


Table  33.  Main  Effect  Test:  System  Size. 


L-M 

Not  L-M 

Row  Totals 

k  =  5 

8 

8 

16 

k=l0 

10 

6 

16 

Column  Totals 

18 

14 

32 

P- Value  =  0.3612  r  =  0.016 


70 


Table  34.  Main  Effect  Test:  Sample  Size. 


L-M 

Not  L-M 

Row  Totals 

Hi  =  4 

9 

7 

16 

n,  =  10 

9 

7 

16 

Column  Totals 

18 

14 

32 

P- Value  =  0.6388  r  =  0.000 


Table  35.  Main  Effect  Test:  Ri  Levels. 


L-M 

Not  L-M 

Row  Totals 

equal 

9 

7 

16 

unequal 

9 

HDH 

16 

Column  Totals 

18 

14 

32 

P- Value  =  0.6388  t  =  0.000 


Table  36.  Main  Effect  Test:  Prior  Information  Strength. 


L-M 

Row  Totals 

uniformly  weak 

10 

6 

16 

8 

8 

16 

Column  Totals 

18 

14 

32 

P- Value  =  0.3612  t  =  0.016 
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Table  37.  Main  Effect  Test;  Prior  Information  Accuracy. 


L-M 

Not  L-M 

Row  Totals 

uniformly  inaccurate 

16 

0 

16 

uniformly  accurate 

2 

14 

16 

Column  Totals 

18 

14 

32 

P- Value  =  0.0000  r  =  0.778 


Table  38.  Interaction  Test;  System  Size  and  Sample  Size. 


L-M 

Not  L-M 

Row  Totals 

4 

4 

8 

k  =  5  and  =  10 

4 

4 

8 

5 

3 

8 

fc  =  10  and  Tii  =  10 

5 

3 

8 

Column  Totals 

18 

14 

32 

=  0.508  P-Value  =  0.9171  t  =  0.016 


Table  39.  Interaction  Test;  System  Size  and  Ri  Levels. 


L-M 

Not  L-M 

Row  Totals 

k  =  5  and  Ri  equal 

4 

4 

8 

k  =  5  and  Ri  unequal 

4 

4 

8 

A;  =  10  and  i?,  equal 

5 

3 

8 

fc  =  10  and  Ri  unequal 

5 

3 

8 

Column  Totals 

18 

14 

32 

=  0.508  P-Value  =  0.9171  t  =  0.016 
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Table  40.  Interaction  Test:  System  Size  and  Prior  Information  Strength. 


L-M 

Not  L-M 

Row  Totals 

k  —  b  and  weak  priors 

4 

4 

8 

k  =  b  and  strong  priors 

4 

4 

8 

A:  =  10  and  weak  priors 

6 

2 

8 

A;  =  10  and  strong  priors 

4 

4 

8 

Column  Totals 

18 

14 

32 

=  1.524  P- Value  =  0.6768  t  =  0.048 


Table  41.  Interaction  Test:  System  Size  and  Prior  Information  Accuracy. 


Not  L-M 

Row  Totals 

A:  =  5  and  inaccurate  priors 

8 

0 

8 

k  =  b  and  accurate  priors 

0 

8 

8 

A:  =  10  and  inaccurate  priors 

8 

0 

8 

A;  =  10  and  accurate  priors 

2 

6 

8 

Column  Totals 

18 

14 

32 

=  25.905  P- Value  =  0.0000  r  =  0.810 


Table  42.  Interaction  Test:  Sample  Size  and  72,  Levels. 


L-M 

Not  L-M 

Row  Totals 

ni  =  4  and  Ri  equal 

5 

3 

8 

Tii  =  4  and  Ri  unequal 

4 

4 

8 

n,  =  10  and  Ri  equal 

4 

4 

8 

n,  =  10  and  Ri  unequal 

5 

3 

8 

Column  Totals 

18 

14 

32 

=  0.508  P- Value  =  0.9171  r  =  0.016 
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Table  43.  Interaction  Test:  Sample  Size  and  Prior  Information  Strength. 


L-M 

Not  L-M 

Row  Totals 

n,  =  4  and  weak  priors 

5 

3 

8 

R,  =  4  and  strong  priors 

4 

4 

8 

n,  =  10  and  weak  priors 

5 

3 

8 

Rj  =  10  and  strong  priors 

4 

4 

8 

Column  Totals 

18 

14 

32 

=  0.508  P- Value  =  0.9171  r  =  0.016 


Table  44.  Interaction  Test:  Sample  Size  and  Prior  Information  Accuracy. 


L-M 

Row  Totals 

n,  =  4  and  inaccurate  priors 

8 

0 

8 

n,  =  4  and  accurate  priors 

1 

7 

8 

Rj  =  10  and  inaccurate  priors 

8 

0 

8 

Ri  =  10  and  accurate  priors 

1 

7 

8 

Column  Totals 

18 

14 

32 

X^  =  24.889  P-Value  =  0.0000  t  =  0.778 


Table  45.  Interaction  Test:  Ri  Levels  and  Prior  Information  Strength. 


L-M 

Hi  equal  and  weak  priors 

5 

3 

8 

Ri  equal  and  strong  priors 

4 

4 

8 

Ri  unequal  and  weak  priors 

5 

3 

8 

Ri  unequal  and  strong  priors 

4 

4 

8 

Column  Totals 

18 

14 

32 

X^  =  0.508  P- Value  =  0.9171  r  =  0.016 


Table  46.  Interaction  Test;  Ri  Levels  and  Prior  Information  Accuracy. 


L-M 

Not  L-M 

Row  Totals 

Ri  equal  and  inaccurate  priors 

8 

0 

8 

Ri  equal  and  accurate  priors 

1 

7 

8 

Ri  unequal  and  inaccurate  priors 

8 

0 

8 

Ri  unequal  and  accurate  priors 

1 

7 

8 

Column  Totals 

18 

14 

32 

=  24.889  P- Value  =  0.0000  r  =  0.778 


Table  47.  Interaction  Test:  Prior  Information  Strength  and  Accuracy. 


L-M 

Not  L-M 

Row  Totals 

weak,  inaccurate  priors 

8 

0 

8 

weak,  accurate  priors 

2 

6 

8 

strong,  inaccurate  priors 

8 

0 

8 

strong,  accurate  priors 

0 

8 

8 

Column  Totals 

18 

14 

32 

=  25.905  P- Value  =  0.0000  r  =  0.810 
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